Dynamic Reliability Analysis of Rotary Steering Drilling System

Vibration and high shock are major factors in the failure of downhole tools. It is important to study the causes of vibration and shock formation to prevent failure of the drillstring and bottom hole assembly (BHA). At present, it is generally recognized that the vibration of drillstring is the main reason for the failure, especially the lateral vibration. In this paper, the bottom tool of Rotary Steering Drilling System (RSS) calculation model was established based on the secondary development of ABAQUS software. Starting from the initial configuration of drilling tool, considering the contact impact of drilling tool and borehole wall, the dynamic excitation of guide mechanism and the drilling pressure, torque, rotational speed, gravity, buoyancy, drilling fluid damping. The dynamic characteristics of the inherent frequency and dynamic stress of the bottom hole assembly (BHA) were calculated and analyzed, and risk assessment method based on the quantitative vibration intensity was established. The reliability of typical drilling tool is evaluated, which provides a reference for the optimization design of BHA of Rotary Steering Drilling System.

1 Introduction

Drillstring and downhole tool failure usually results from failing to control one or more of the vibration mechanisms. The drillstring motion state is affected by several factors and interaction among them, like drilling load, drillstring itself, drilling fluid, well wall and bottom rock (Wang et al., 2015). The vibration of the drillstring can be divided into three primary modes: axial, lateral (also referred to as transverse, or bending), and torsional modes, with different destructive natures (Ritto et al., 2013; Liao et al., 2012). In the axial mode, the vibration is a longitudinal movement along the drillstring. In the transverse mode, the vibration is a lateral movement that causes bending or bending of the component, resulting in a stress reversal in which one side is in a different tension state than the other side. In the torsional mode, the vibration is resistance to the rotation resulting in twisting as torque is applied to overcome resistance (Nave and Dvorkin, 2015; Ghasemloonia et al., 2014; Ritto et al., 2013).

Wilson and Heisig (2015) found that lateral vibration is more severe than vertical vibration, because the bending stress is high while vertical stress caused by axial vibration is not severe. Lateral vibration is very easy to happen in side tracking, and the lowest resonant frequency depends on properties of drillstring material and drilling fluid. Qu (1995) discussed lateral vibration under the effect of drilling fluid inside and outside by using solid-liquid coupling vibration theory, he concludes that drilling fluid can be equivalent to generalized distributive mass, and the effect of drilling fluid cannot be neglected. Shor et al. (2015) utilizing finite element method get conclusion that mass unbalance has a great influence on lateral vibration, and when the drillstring vibrates at the inherent frequency, the friction between the drillstring and well wall becomes more severe. Drillstring dynamics study are mainly by method of theoretical computation due to the difficulty of experiments. Millheim et al. (1978) lay the foundation of using finite element in dynamics analysis. The finite element computer program to analyze the behavior of bottom-hole assemblies of drill collars, stabilizers and drill bit is discussed. Wilson and Heisig (2015), Schmalhorst and Neubert (2003), Jogi et al. (2002) and Yigit et al. (1996) utilized both differential equation solution and finite element method to analysis the drillstring dynamics, and some of their work has formed commercial software. In most cases analytical models are used to get a fast overview of natural axial and torsional frequencies of drillstrings in straight boreholes. To develop practical and commercial computing software about drillstring dynamics, the model established must take lithology, bit type, drill assembly, damping and contact into consideration, meanwhile through a lot of experiments and experiment data analysis, comparison with theoretical computing results, various parameters or boundary conditions in theoretical calculation can be corrected. Therefore, it takes a lot of hard work to develop software for drillstring dynamics. Based on this situation, secondary development based on commercial finite element software is one feasible technical approach, which reduces large amount of work for calculation method and program verification.

Several studies have been devoted to better understand the full dynamics of the rotary drilling systems (Yigit and Christoforou, 1996, 1998; Khulief and Al-Naser, 2005; Sampaio et al., 2007; Kapitaniak et al., 2015). Gupta and Wahi (2016) developed a global dynamic of the coupled axial-torsional vibrations incorporating the possibility of bit-bounce. Huang et al. (2015) established a generalized quasi-static model of drillstring system, the forward model and inversion model are further combined into a unity. Lian et al. (2015) established a nonlinear dynamics model to investigate the motion behavior of drillstring in gas drilling of horizontal wells, a finite element model was also established and the buckling and contact of drillstring was analyzed. However, none of these studies have reported the dynamical motion state of the bottom hole assembly (BHA) in Rotary Steerable Drilling System (RSS) (Jones et al., 2016), the whole drillstring is rotated from the surface by a hydraulically driven top drive. The RSS has push-the-bit tools (Warren, 1998) which using pads on the outside of the tool which press against the well bore thereby causing the bit to press on the opposite side causing a direction change. The pads of the implementing agencies in RSS constantly pushed against the borehole wall, making bottom hole by a cycle of nonlinear damping force, which is lead to the bottom drilling tool movement of chaos and disorder (Wang et al., 2014; Xue et al., 2015, 2019).

Statistics study (Lin et al., 2016) shows that drillstring failure and cross-over sub failure account for approximately 54 % of drilling tool failures. Among the drillstring thread fracture incidents, 70 % of them occurred at box thread and 30 % occurred at pin thread. The causes for drillstring failure encountered during drilling practice have been systematically analyzed based on the statistical data (Wang et al., 2011; Abdo, 2006; Abdo and Abouelsoud, 2011). Fault tree of drillstring failure has been subsequently constructed (Wang et al., 2011).

In this paper, the bottom tool of RSS calculation model was established based on the secondary development of ABAQUS software. Starting from the initial configuration of drilling tool, considering the contact impact of drilling tool and borehole wall, the dynamic excitation of guide mechanism and the drilling pressure, torque, rotational speed, gravity, buoyancy, drilling fluid damping. The dynamic characteristics of the inherent frequency and dynamic stress of the BHA were calculated and analyzed, and risk assessment method based on the quantitative vibration intensity was established. The reliability of typical drilling tool is evaluated, which provides a reference for the optimization design of drilling tools.

2 Mathematical Model of Drillstring Dynamics

Drillstring actual working condition is very complex, and it's hard to simulate and analyze the actual working condition precisely. Due to the complexity of the drilling operation, the mathematical calculation model cannot effectively explain the fatigue failure of the drillstring caused by the cyclic stress due to vibration and high shock. Since the dynamics of the drillstring involve wide vibration profiles including axial, lateral and torsional modes, the mathematical modeling of such long rotating components is highly nonlinearity.

Therefore, secondary development based on commercial finite element software is a feasible technical solution. In this paper, the drillstring be simplified by basic assumptions when establish a vibration analysis model in ABAQUS software. The BHA calculation model of the RSS is shown in Fig. 1, making the following assumptions in the calculation process of the entire dynamic system.

• The borehole section is circular and the well wall is rigid;

• The axis of drillstring is at the same line with well axis before drillstring moves, and the drillstring node has three-dimensional freedom;

• The bottom drillstring is considered as a beam, whose one end is hinging support and the other end is movable hinge support. The drill bit and the stabilizer are movable hinge support, push contact is one-way movable support;

• The drillstring is considered as homogeneous, small deformation and elastic beam, the joint is neglected, and drillstring displacement is restrained by wellbore;

• The drilling fluid influence on drillstring vibration characteristics is considered.

Figure 1 Dynamics model of Rotary Steering Drilling System (RSS).

Dynamic equilibrium equation of drilling dynamics system can be expressed as:

$$\begin{array}{}\text{(1)}& {\displaystyle }\mathbf{M}\ddot{U}+\mathbf{C}\dot{U}+\mathbf{K}U-P=\mathrm{0}\end{array}$$

wherein U is generalized displacement, $\dot{U}$ is generalized velocity, $\ddot{U}$ is generalized acceleration, P is external force, M is mass matrix, K is stiffness matrix, C is damping matrix.

Considering that the inside of the drill is filled with liquid and the outer annulus also has drilling fluid, the unit of drillstring could have the ability to withstand greater stress and deflection. The stiffness matrix of the unit is shown in Eq. (2) (the equation details can be found in the Supplement).

$$\begin{array}{ll}\text{(2)}& {\displaystyle }& {\displaystyle }{\mathbf{K}}_e={\displaystyle }& {\displaystyle }\left[\begin{array}{cccccccccccc}AE/L& & & T_{\mathrm{T}}& & & AE/L& & & -T_{\mathrm{T}}& & \\ & a_Z& & & & c_Z& & -a_Z& & & & c_Z\\ & & a_Y& & -c_Y& & & & -a_Y& & -c_Y& \\ T_{\mathrm{T}}& & & GJ/L& & & -T_{\mathrm{T}}& & & -GJ/L& & \\ & & -c_Y& & e_Y& & & & c_Y& & f_Y& \\ & c_Z& & & & e_Z& & -c_Z& & & & f_Z\\ AE/L& & & -T_{\mathrm{T}}& & & AE/L& & & T_{\mathrm{T}}& & \\ & -a_Z& & & & -c_Z& & a_Z& & & & -c_Z\\ & & -a_Y& & c_Y& & & & a_Y& & c_Y& \\ -T_{\mathrm{T}}& & & -GJ/L& & & T_{\mathrm{T}}& & & GJ/L& & \\ & & -c_Y& & f_Y& & & & c_Y& & e_Y& \\ & c_Z& & & & f_Z& & -c_Z& & & & e_Z\end{array}\right]\end{array}$$

wherein, A is cross-sectional area, E is Young's modulus, L is unit length, G is shear modulus, J is torsional moment of inertia, if I_x=0, then define J_x, if I_x≠0, then define I_i, where I_x is input torsional moment of inertia about the x axis, I_y is moment of inertia about y axis, I_z is moment of inertia about z axis; $J_x=I_y+I_z$ is polar moment of inertia. Then we can define $a_Z=a(I_Z,{\mathit{\varphi }}_Y)$, $a_Y=a(I_Y,{\mathit{\varphi }}_Z)$, $b_Z=b(I_Z,{\mathit{\varphi }}_Y)$, $c_z=c(I_z,{\mathit{\phi }}_y)$, …, $f_Z=f(I_Z,{\mathit{\varphi }}_Y)$, $f_Y=f(I_Y,{\mathit{\varphi }}_Z)$, the following equations can be obtained.

$$\begin{array}{}\text{(3)}& {\displaystyle }& {\displaystyle }a(I,\mathrm{\Phi })={\displaystyle \frac{\mathrm{12}EI}{L^{\mathrm{3}}(\mathrm{1}+\mathrm{\Phi })}}\text{(4)}& {\displaystyle }& {\displaystyle }c(I,\mathrm{\Phi })={\displaystyle \frac{\mathrm{6}EI}{L^{\mathrm{2}}(\mathrm{1}+\mathrm{\Phi })}}\text{(5)}& {\displaystyle }& {\displaystyle }e(I,\mathrm{\Phi })={\displaystyle \frac{(\mathrm{4}+\mathrm{\Phi })EI}{L(\mathrm{1}+\mathrm{\Phi })}}\text{(6)}& {\displaystyle }& {\displaystyle }f(I,\mathrm{\Phi })={\displaystyle \frac{(\mathrm{2}-\mathrm{\Phi })EI}{L(\mathrm{1}+\mathrm{\Phi })}}\text{(7)}& {\displaystyle }& {\displaystyle }{\mathrm{\Phi }}_Y={\displaystyle \frac{\mathrm{12}E{I}_Z}{G{A}_Z^{\mathrm{S}}{L}^{\mathrm{2}}}}\text{(8)}& {\displaystyle }& {\displaystyle }{\mathrm{\Phi }}_Z={\displaystyle \frac{\mathrm{12}E{I}_Y}{G{A}_Y^{\mathrm{S}}{L}^{\mathrm{2}}}}\end{array}$$

wherein I_i is the normal moment of inertia in the i direction, $A_i^{\mathrm{S}}$ is the shear area perpendicular to the direction i, where $i=A/F_i^{\mathrm{S}}$, $F_i^{\mathrm{S}}$ is the shear coefficient. Then we can get the torque T_T as following:

$$\begin{array}{}\text{(9)}& {\displaystyle }{T}_{\mathrm{T}}=\left(\begin{array}{c}\mathrm{0}\\ \\ {\displaystyle \frac{G_{\mathrm{T}}(D_{\mathrm{o}}^{\mathrm{3}}-D_{\mathrm{i}}^{\mathrm{3}})}{L}}\end{array}\right.\end{array}$$

wherein, the standard option for torque is 0, and the twisted pair for torque imbalance is $G_{\mathrm{T}}(D_{\mathrm{o}}^{\mathrm{3}}-D_{\mathrm{i}}^{\mathrm{3}})/L$; G_T is the twisted pair tensile stiffness constant; D_i is the in-tube diameter, its value is equal to D_o−2t_w; t_w is the wall thickness, mm; D_o is the outer diameter, mm.

The mass matrix of the unit is shown in Eq. (10) (the equation details can be found in the Supplement), and its mass matrix is similar to the mass matrix of ordinary beam elements. Only some elements need to be revised, multiplied by the coefficient M_a∕M_t.

$$\begin{array}{ll}\text{(10)}& {\displaystyle }& {\displaystyle }{\mathbf{M}}_e={\displaystyle }& {\displaystyle }{M}_t\left[\begin{array}{cccccccccccc}{M}_a/\mathrm{3}{M}_t& & & & & & E{M}_{a}A/M_{t}L& & & & & \\ & A_Z& & & & C_Z& & B_Z& & & & -D_Z\\ & & A_Y& & -C_Y& & & & B_Y& & D_Y& \\ & & & J_X{M}_a/\mathrm{3}{M}_{t}A& & & & & & -G{M}_{a}J/M_{t}L& & \\ & & -C_Y& & E_Y& & & & -D_Y& & F_Y& \\ & C_Z& & & & E_Z& & D_Z& & & & F_Z\\ M_a/\mathrm{6}{M}_t& & & & & & M_a/\mathrm{3}{M}_t& & & & & \\ & B_Z& & & & D_Z& & A_Z& & & & C_Z\\ & & B_Y& & -D_Y& & & & A_Y& & -C_Y& \\ & & & J_X{M}_a/\mathrm{6}{M}_{t}A& & & & & & J_X{M}_a/\mathrm{3}{M}_{t}A& & \\ & & D_Y& & F_Y& & & & -C_Y& & E_Y& \\ & -D_Z& & & & F_Z& & C_Z& & & & E_Z\end{array}\right]\end{array}$$

wherein,

$$\begin{array}{}\text{(11)}& {\displaystyle }& {\displaystyle }{M}_a=(m_{\mathrm{w}}+m_{\mathrm{int}})L\text{(12)}& {\displaystyle }& {\displaystyle }{M}_t=(m_{\mathrm{w}}+m_{\mathrm{int}}+m_{\mathrm{ins}}+m_{\mathrm{add}})L\text{(13)}& {\displaystyle }& {\displaystyle }{m}_{\mathrm{w}}=\mathit{\rho }{\displaystyle \frac{\mathit{\pi }(D_{\mathrm{o}}^{\mathrm{2}}-D_{\mathrm{i}}^{\mathrm{2}})}{\mathrm{4}}}\text{(14)}& {\displaystyle }& {\displaystyle }{m}_{\mathrm{add}}=C_M{\mathit{\rho }}_{\mathrm{f}}{\displaystyle \frac{\mathit{\pi }{D}_{\mathrm{o}}^{\mathrm{2}}}{\mathrm{4}}}\end{array}$$

wherein, m_w is the mass of the drillstring per unit length; m_int is the mass of the internal fluid per unit length; m_ins is the quality of the external additional layer, set m_ins=0; m_add is the additional mass of external fluid; ρ is the density of drillstring; C_M is the additional mass factor for external fluid; ρ_f is fluid density.

Then we can define $A_Z=A(r_Z,{\mathit{\varphi }}_Y)$, $A_Y=A(r_Y,{\mathit{\varphi }}_Z)$, $B_Z=B(r_Z,{\mathit{\varphi }}_Y)$, … $F_Z=F(r_Z,{\mathit{\varphi }}_Y)$, $F_Y=F(r_Y,{\mathit{\varphi }}_Z)$, radius of rotation $r_y=\sqrt{I_{yy}/A}$, and $r_z=\sqrt{I_{zz}/A}$. Detailed calculation expression is as follows:

$$\begin{array}{}\text{(15)}& {\displaystyle }& {\displaystyle }A(r,\mathit{\phi })={\displaystyle \frac{\frac{\mathrm{13}}{\mathrm{35}}+\frac{\mathrm{7}}{\mathrm{10}}\mathit{\phi }+\frac{\mathrm{1}}{\mathrm{3}}{\mathit{\phi }}^{\mathrm{2}}+\frac{\mathrm{6}}{\mathrm{5}}(r/L{)}^{\mathrm{2}}}{(\mathrm{1}+\mathit{\phi }{)}^{\mathrm{2}}}}\text{(16)}& {\displaystyle }& {\displaystyle }B(r,\mathit{\phi })={\displaystyle \frac{\frac{\mathrm{9}}{\mathrm{70}}+\frac{\mathrm{3}}{\mathrm{10}}\mathit{\phi }+\frac{\mathrm{1}}{\mathrm{6}}{\mathit{\phi }}^{\mathrm{2}}+\frac{\mathrm{6}}{\mathrm{5}}(r/L{)}^{\mathrm{2}}}{(\mathrm{1}+\mathit{\phi }{)}^{\mathrm{2}}}}\text{(17)}& {\displaystyle }& {\displaystyle }C(r,\mathit{\phi })={\displaystyle \frac{\left[\frac{\mathrm{11}}{\mathrm{210}}+\frac{\mathrm{11}}{\mathrm{120}}\mathit{\phi }+\frac{\mathrm{1}}{\mathrm{24}}{\mathit{\phi }}^{\mathrm{2}}+(\frac{\mathrm{1}}{\mathrm{10}}-\frac{\mathrm{1}}{\mathrm{2}}\mathit{\phi })(r/L{)}^{\mathrm{2}}\right]L}{(\mathrm{1}+\mathit{\phi }{)}^{\mathrm{2}}}}\text{(18)}& {\displaystyle }& {\displaystyle }D(r,\mathit{\phi })={\displaystyle \frac{\left[\frac{\mathrm{13}}{\mathrm{420}}+\frac{\mathrm{3}}{\mathrm{40}}\mathit{\phi }+\frac{\mathrm{1}}{\mathrm{24}}{\mathit{\phi }}^{\mathrm{2}}+(\frac{\mathrm{1}}{\mathrm{10}}-\frac{\mathrm{1}}{\mathrm{2}}\mathit{\phi })(r/L{)}^{\mathrm{2}}\right]L}{(\mathrm{1}+\mathit{\phi }{)}^{\mathrm{2}}}}\text{(19)}& {\displaystyle }& {\displaystyle }E(r,\mathit{\phi })={\displaystyle \frac{\left[\frac{\mathrm{1}}{\mathrm{105}}+\frac{\mathrm{1}}{\mathrm{60}}\mathit{\phi }+\frac{\mathrm{1}}{\mathrm{120}}{\mathit{\phi }}^{\mathrm{2}}+(\frac{\mathrm{2}}{\mathrm{15}}+\frac{\mathrm{1}}{\mathrm{6}}\mathit{\phi }+\frac{\mathrm{1}}{\mathrm{3}}{\mathit{\phi }}^{\mathrm{2}})(r/L{)}^{\mathrm{2}}\right]L^{\mathrm{2}}}{(\mathrm{1}+\mathit{\phi }{)}^{\mathrm{2}}}}\text{(20)}& {\displaystyle }& {\displaystyle }F(r,\mathit{\phi })={\displaystyle \frac{\left[\frac{\mathrm{1}}{\mathrm{140}}+\frac{\mathrm{1}}{\mathrm{60}}\mathit{\phi }+\frac{\mathrm{1}}{\mathrm{120}}{\mathit{\phi }}^{\mathrm{2}}+(\frac{\mathrm{1}}{\mathrm{30}}+\frac{\mathrm{1}}{\mathrm{6}}\mathit{\phi }-\frac{\mathrm{1}}{\mathrm{6}}{\mathit{\phi }}^{\mathrm{2}})(r/L{)}^{\mathrm{2}}\right]L^{\mathrm{2}}}{(\mathrm{1}+\mathit{\phi }{)}^{\mathrm{2}}}}\end{array}$$

The natural frequency (Jogi et al., 2002) can be obtained by analyzing the dynamic response of the structure when the external load is zero (in the dynamic equilibrium equation P=0). In Abaqus software, the vibration mode is calculated for the undamped system. In this case, the equation of motion can be expressed as:

$$\begin{array}{}\text{(21)}& {\displaystyle }\mathbf{K}-{\mathit{\omega }}^{\mathrm{2}}\mathbf{M}=\mathrm{0}\end{array}$$

wherein ω is the natural frequency of the drillstring system.

The n-order natural frequency and the n-order natural mode of the drillstring can be obtained by the Eq. (21), and calculation of the natural frequency and the natural mode is mathematically a problem of calculate all the eigenvalues of the matrix. For the vibration damage caused by the general engineering structure, it usually only occurs in the lower frequency range, so only a few lower order natural frequencies will be obtained.

However, most engineering problems still involve damping, although the damping may be small. The relationship of natural frequency between the damped and the undamped is:

$$\begin{array}{}\text{(22)}& {\displaystyle }{\mathit{\omega }}_{\mathrm{d}}=\mathit{\omega }\sqrt{\mathrm{1}-{\mathit{\xi }}^{\mathrm{2}}}\end{array}$$

wherein, ω_d is the natural frequency of the system when damping, ξ is the damping ratio.

Based on the above analysis and Hamilton's principle, the dynamic model of the drillstring system can be derived as follows:

$$\begin{array}{}\text{(23)}& {\displaystyle }\mathbf{M}\ddot{u}+F_{\mathrm{F}}(u,\dot{u})+F_{\mathrm{w}}(u,\dot{u})+F_{\mathrm{G}}\left(u\right)=R+F_{\mathrm{E}}(u,\dot{u},t)\end{array}$$

Where, u is node displacement or rotation, rad; M is the mass matrix; F_F is drilling fluid distributive force, kN m⁻¹; F_w is well bore contact force, kN; F_G is nonlinear elastic force, kN; R is static force (gravity, buoyancy, WOB, etc.), kN; F_E is exciting force (uneven quality force, axial bit load, lateral bit load), kN. The load form and load condition need to be selected and specialized. Due to working characteristics of RSS, all external forces, such as push force, well wall contact force, drilling fluid distributive force, nonlinear elastic force, are taken into consideration besides bit load, gravity and buoyancy.

3 Reliability Analysis Method

3.1 Risk assessment method based on quantitative vibration intensity

The drillstring failure and cross-over sub failure account for approximately 54 % of drilling tool failures. Drillstring and downhole tool failure usually results from failing to control one or more of the vibration mechanisms. We established a risk assessment method based on quantitative vibration intensity, the evaluation process is shown in Fig. 2.

Figure 2 Analysis process of drillstring vibration intensity.

As show in Fig. 2, according to the initial designed speed and drilling tool combination, the calculation model is established. Firstly, the critical speed is calculated, and the range of relevant frequencies is taken to determine the danger control point. The drilling tool is not allowed to be attached to the critical speed. For a specific speed, it may be between the two critical speeds, the low critical speeds are expressed as n_c1, and the high critical speeds are expressed as n_c2, if the design speed is satisfied $\mathrm{1.3}{n}_{\mathrm{c}\mathrm{1}}\le n\le \mathrm{0.75}{n}_{\mathrm{c}\mathrm{2}}$, no resonance will occur and the point control requirements will be satisfied. If the above conditions are not satisfied, the recommended speed value is given as the intermediate value of the two-stage critical speed, that is $n^{\prime }=\mathrm{0.5}(n_{\mathrm{c}\mathrm{1}}+n_{\mathrm{c}\mathrm{2}})$.

From the calculation results above, the length of drillstring is large, the natural frequency is small, and sometimes the recommended value is difficult to stay away from the resonance speed. Of course, because of the existence of the wellbore, even if the drillstring speed is close to the critical speed, there will be no infinite amplitude, but the collision between the drillstring and the wellbore will be intensified. When the rotating speed is close to the natural frequency, under the action of periodic excitation load, the drillstring is usually in some violent vibration state, and the drillstring has a certain dynamic stress level. Further stress level judgment is needed to determine the reliability of the recommended speed.

The mechanical properties of the drillstring need to be determined when determine stress intensity, that include the tensile strength σ_b, yield limit σ_s, symmetrical cyclic fatigue limit σ₋₁ and allowable stress amplitude σ_a which considering the effects of mean stress. According to the dynamic model, the vibration stress will be obtained, and the time series of the equivalent stress will be obtained, so the transient stress σ, the equivalent stress amplitude σ_da and the equivalent average stress σ_dn can be obtained. And by contrasting allowable stress amplitude σ_a, the fatigue limit σ₋₁, and yield limit σ_s, the condition is defined as “safety” when equivalent stress amplitude σ_da is low, less than the revised allowable fatigue stress amplitude σ_a; The condition is defined as “warning” when equivalent stress amplitude σ_da exceeds revised allowable fatigue stress amplitude σ_a, but less than the fatigue limit σ₋₁, and should be handled carefully, it is not allowed to occur that the equivalent stress amplitude exceeds the fatigue limit σ₋₁ or the instantaneous stress σ exceeds the yield strength σ_s. Therefore, it is necessary for the construction personnel to make quick treatment and control the vibration intensity.

3.2 Calculation method and judgment principle

3.2.1 Calculation of equivalent stress amplitude and equivalent average stress

When the fatigue strength check with known three-way stress, the multi-direction stress is often converted into equivalent stress in engineering, and then the fatigue check is carried out according to the stress state. According to a large number of experiments, the deformation energy intensity theory combines the multidirectional stress state with the unidirectional stress state, which is more in line with the practical theory. According to the deformation energy intensity theory, the multidirectional stress will transform into a unidirectional force.

Then the equivalent stress amplitude σ_da can be defined:

$$\begin{array}{}\text{(24)}& {\displaystyle }{\mathit{\sigma }}_{\mathrm{da}}=\sqrt{{\displaystyle \frac{{\left({\mathit{\sigma }}_{\mathrm{a}\mathrm{1}}-{\mathit{\sigma }}_{\mathrm{a}\mathrm{2}}\right)}^{\mathrm{2}}+{\left({\mathit{\sigma }}_{\mathrm{a}\mathrm{2}}-{\mathit{\sigma }}_{\mathrm{a}\mathrm{3}}\right)}^{\mathrm{2}}+{\left({\mathit{\sigma }}_{\mathrm{a}\mathrm{3}}-{\mathit{\sigma }}_{\mathrm{a}\mathrm{1}}\right)}^{\mathrm{2}}}{\mathrm{2}}}}\end{array}$$

The equivalent mean stress σ_dm will be defined as:

$$\begin{array}{}\text{(25)}& {\displaystyle }{\mathit{\sigma }}_{\mathrm{dm}}=\sqrt{{\displaystyle \frac{{\left({\mathit{\sigma }}_{\mathrm{m}\mathrm{1}}-{\mathit{\sigma }}_{\mathrm{m}\mathrm{2}}\right)}^{\mathrm{2}}+{\left({\mathit{\sigma }}_{\mathrm{m}\mathrm{2}}-{\mathit{\sigma }}_{\mathrm{m}\mathrm{3}}\right)}^{\mathrm{2}}+{\left({\mathit{\sigma }}_{\mathrm{m}\mathrm{3}}-{\mathit{\sigma }}_{\mathrm{m}\mathrm{1}}\right)}^{\mathrm{2}}}{\mathrm{2}}}}\end{array}$$

wherein σ_a1, σ_a2 and σ_a3 are the main stress amplitude, MPa; σ_m1, σ_m2 and σ_m3 are the average magnitude of principal stress, MPa.

If the time series of equivalent stress is obtained by measurement, its equivalent stress σ_da amplitude and equivalent average stress σ_dn can be obtained.

where equivalent stress amplitude expressed as:

$$\begin{array}{}\text{(26)}& {\displaystyle }{\mathit{\sigma }}_{\mathrm{da}}={\displaystyle \frac{{\mathit{\sigma }}_{max}-{\mathit{\sigma }}_{min}}{\mathrm{2}}}\end{array}$$

Equivalent mean stress expressed as:

$$\begin{array}{}\text{(27)}& {\displaystyle }{\mathit{\sigma }}_{\mathrm{dm}}={\displaystyle \frac{{\mathit{\sigma }}_{max}+{\mathit{\sigma }}_{min}}{\mathrm{2}}}\end{array}$$

3.2.2 Calculation of mechanical performance parameters

To simplified the calculation model, combined with the fatigue endurance limit test of drillstring, the approximate calculation formula of the cycle endurance limit of drillstring under the condition of underground corrosion was adopted. The fatigue limit σ₋₁ under symmetrical circulation is shown in Eq. (28).

$$\begin{array}{}\text{(28)}& {\displaystyle }{\mathit{\sigma }}_{-\mathrm{1}}=\mathrm{0.3}\left({\mathit{\sigma }}_{\mathrm{s}}+{\mathit{\sigma }}_{\mathrm{b}}\right)\end{array}$$

wherein σ_s is the yield strength, MPa, σ_b is tensile limit, MPa.

The fatigue allowable stress amplitude formula for under symmetrical circulation expressed as:

$$\begin{array}{}\text{(29)}& {\displaystyle }{\mathit{\sigma }}_{-\mathrm{1}}={\displaystyle \frac{\mathit{\beta }{\mathit{\epsilon }}_{\mathit{\sigma }}}{K_{\mathit{\sigma }}S}}{\mathit{\sigma }}_{-\mathrm{1}}\end{array}$$

wherein, σ₋₁ is the allowable stress amplitude, MPa; β is surface mass coefficient; K_σ is effective stress concentration factor; ε_σ is size coefficient; σ₋₁ is fatigue limit under symmetric circulation, MPa; S is the safety factor. According to the actual situation of drillstring, we set K_σ=1, ε_σ=0.74, β=0.9, S=2.

It is necessary to modify the allowable fatigue stress value when considering the effect of average stress according to the previous calculation. Goodman line is often used for modification in engineering, which is relatively simple and the calculation results tend to be safe, so it is widely used (Budynas and Nisbett, 2014). The Goodman's line calculation formula can be expressed as:

$$\begin{array}{}\text{(30)}& {\displaystyle }{\mathit{\sigma }}_{\mathrm{a}}={\mathit{\sigma }}_{\mathrm{D}}\left(\mathrm{1}-{\displaystyle \frac{{\mathit{\sigma }}_{\mathrm{m}}}{{\mathit{\sigma }}_{\mathrm{b}}}}\right)={\mathit{\sigma }}_{-\mathrm{1}}\left(\mathrm{1}-{\displaystyle \frac{{\mathit{\sigma }}_{\mathrm{dm}}}{{\mathit{\sigma }}_{\mathrm{b}}}}\right)\end{array}$$

wherein, σ_a is allowable stress amplitude considering the effects of mean stress, σ_D is the allowable stress amplitude under symmetrical circulation, equal to σ₋₁, MPa; σ_m is equivalent mean stress, MPa; σ_b is tensile strength, MPa.

3.2.3 Principle of judgment

The judgment principle of vibration stress intensity is defined in this paper as following:

1. σ_da≤σ_a, drillstring vibration intensity is small, can be identified as safely.

2. ${\mathit{\sigma }}_{\mathrm{a}}\le {\mathit{\sigma }}_{\mathrm{da}}\le {\mathit{\sigma }}_{-\mathrm{1}}$, the vibration intensity of drillstring is large, can be identified as warning.

3. σ≥σ_s or ${\mathit{\sigma }}_{\mathrm{da}}\le {\mathit{\sigma }}_{-\mathrm{1}}$, the vibration intensity of drillstring is too high, can be identified as dangerous.

3.3 Secondary development software

In this paper, we established calculation program of the RSS based on the secondary development of Abaqus software. The modeling starts from the drill bit, input the length, outer diameter and inner diameter of each section in sequence according to the different cross-sections of drilling tool combination. In the column of material properties, the elastic modulus, Poisson's ratio, and density of materials are respectively input, where the unit of elastic modulus is MPa and the unit of density is t mm⁻³.

Load column input weight on bit, rotary torque, slap backup force and impact force, according to the actual situation input weight on bit, rotary table torque, slap pushing and collision force, the naming principles of section number is the number under push position. At the same time input the drillstring rotary speed in the working process, the unit is r min⁻¹, and then meshing, at last in the task volume input file name and the number of CPU, and then click OK, begin to automatically calculate the dynamic stress of drillstring assembly. In the column of material mechanics parameters, input material yield strength and tensile strength, in the column of calculation results, input maximum instantaneous stress, the maximum and the minimum fluctuation stress, and then click OK, appear evaluation results of assembly dynamic stress.

Calculation process in the software is shown as in Fig. 3. Firstly, by using Abaqus model input, establish drillstring dynamics model, set analysis steps, define contact and impose static load, and simulate quasi-static deformation process of the contract model between drillstring and wellbore which is under axial pressure, and obtain drillstring initial configuration. Analyze the buckling results, the initial displacement is obtained and the deformation model is imported for further dynamic analysis. On the one hand, the collision calculation is carried out to obtain the collision contact force function. On the other hand, the natural frequency of the model is calculated to determine the damping coefficient of the system. The calculated collision contact force function and damping coefficient are introduced into the model, and other load functions such as speed, torque, gravity, buoyancy and pushing force are applied. Submit the work for dynamic calculation. During the calculation, the boundary changes are monitored in real time, and the boundary corrections are made in case of collision. The time history of dynamic stress under various excitation loads and the deformation, velocity and acceleration of drillstring are obtained.

Figure 3 Modal analysis calculation process.

4 Results

4.1 Initial conditions for calculation

BHA design parameters: Drill bit (φ215.9 mm), rotary steerable tool (φ197 mm), drillstring (φ127 mm), undersized centralizer (φ210 mm), flexible sub (φ172 mm), drill collar (φ172 mm, length 6 m), non-magnetic drill pipe (φ127).

The elastic modulus of drill assembly is 206 GPa, Poisson's ratio is 0.3, drillstring material density is 7850 kg m⁻³, and drilling fluid density is 1150 kg m⁻³. Well bore diameter is 216 mm. The rotary table torque applied to the drillstring is 10 kN m, weight on bit is 80 kN, and excitation load is applied to the pushing block. The whole drillstring is affected by gravity, buoyancy and damping force of drilling fluid.

When different angular velocity is applied to the drillstring, the corresponding period will be changed. The drill assembly has one end hinge branch, the other end movable hinge, the centralizer movable hinge, the shaft wall solid branch.

4.2 Dynamic Stress Calculation

The principal stresses in three directions are shown in Fig. 4. The principal stress time series curve satisfies the condition of periodic cycle, but there is error in symmetry.

Figure 4 Time series curve of typical principal stress.

For the first principal stress, the maximum stress ${\mathit{\sigma }}_{\mathrm{1}max}=\mathrm{52.8}$ (MPa), and the minimum stress ${\mathit{\sigma }}_{\mathrm{1}min}=-\mathrm{9.06}$ (MPa), so the main stress amplitude will be calculated as follows.

$${\displaystyle }{\mathit{\sigma }}_{\mathrm{1}\mathrm{a}}={\displaystyle \frac{{\mathit{\sigma }}_{\mathrm{1}max}-{\mathit{\sigma }}_{\mathrm{1}min}}{\mathrm{2}}}={\displaystyle \frac{\mathrm{52.8}-(-\mathrm{9.06})}{\mathrm{2}}}=\mathrm{30.93}\left(\text{MPa}\right)$$

The average stress in the direction of the principal stress amplitude can be calculated:

$${\displaystyle }{\mathit{\sigma }}_{\mathrm{1}\mathrm{m}}={\displaystyle \frac{{\mathit{\sigma }}_{\mathrm{1}max}+{\mathit{\sigma }}_{\mathrm{1}min}}{\mathrm{2}}}={\displaystyle \frac{\mathrm{52.8}+(-\mathrm{9.06})}{\mathrm{2}}}=\mathrm{21.87}\left(\text{MPa}\right)$$

Similarly, the same formula is used to calculate the second and third principal stress, as shown in Table 1. The equivalent stress amplitude σ_da is calculated:

$$\begin{array}{ll}{\displaystyle }{\mathit{\sigma }}_{\mathrm{da}}& {\displaystyle }=\sqrt{{\displaystyle \frac{{\left(\mathrm{30.93}-\mathrm{21.39}\right)}^{\mathrm{2}}+{\left(\mathrm{21.39}-\mathrm{33.59}\right)}^{\mathrm{2}}+{\left(\mathrm{33.59}-\mathrm{30.93}\right)}^{\mathrm{2}}}{\mathrm{2}}}}\\ {\displaystyle }& {\displaystyle }=\mathrm{15.7}\left(\text{MPa}\right)\end{array}$$

equivalent mean stress calculated:

$$\begin{array}{ll}{\displaystyle }{\mathit{\sigma }}_{\mathrm{dm}}& {\displaystyle }=\sqrt{{\displaystyle \frac{\begin{array}{c}{\left(\mathrm{21.87}-(-\mathrm{1.05})\right)}^{\mathrm{2}}+{\left((-\mathrm{1.05})-(-\mathrm{25.52})\right)}^{\mathrm{2}}\\ +{\left((-\mathrm{25.52})-\mathrm{21.87}\right)}^{\mathrm{2}}\end{array}}{\mathrm{2}}}}\\ {\displaystyle }& {\displaystyle }=\mathrm{58.05}\left(\text{MPa}\right)\end{array}$$

Relative to drillstring materials, σ_b=827 (MPa), σ_s=758 (MPa), other coefficients are selected in accordance with Sect. 4.1. We can obtain,

$${\displaystyle }{\mathit{\sigma }}_{-\mathrm{1}}=\mathrm{0.3}\left({\mathit{\sigma }}_{\mathrm{s}}+{\mathit{\sigma }}_{\mathrm{b}}\right)=\mathrm{483.6}\left(\text{MPa}\right)$$

wherein, σ_s is yield strength, MPa; σ_b is tensile strength, MPa.

Then the allowable stress under symmetric circulation is,

$${\displaystyle }{\mathit{\sigma }}_{-\mathrm{1}}={\displaystyle \frac{\mathrm{0.9}\times \mathrm{0.74}}{\mathrm{1.0}\times \mathrm{2}}}\times \mathrm{483.6}=\mathrm{161.0}\left(\text{MPa}\right)$$

According to the above calculation data, the allowable stress amplitude σ_a considering the influence of average stress is,

$$\begin{array}{ll}{\displaystyle }{\mathit{\sigma }}_{\mathrm{a}}& {\displaystyle }={\mathit{\sigma }}_{\mathrm{D}}\left(\mathrm{1}-{\displaystyle \frac{{\mathit{\sigma }}_{\mathrm{m}}}{{\mathit{\sigma }}_{\mathrm{b}}}}\right)={\mathit{\sigma }}_{-\mathrm{1}}\left(\mathrm{1}-{\displaystyle \frac{{\mathit{\sigma }}_{\mathrm{m}}}{{\mathit{\sigma }}_{\mathrm{b}}}}\right)=\mathrm{161}\times \left(\mathrm{1}-{\displaystyle \frac{\mathrm{58.05}}{\mathrm{827}}}\right)\\ {\displaystyle }& {\displaystyle }=\mathrm{149.7}\left(\text{MPa}\right)\end{array}$$

Then, σ is obviously less than σ_s, and σ_da<σ_a, therefore, meets the fatigue strength requirement.

Table 1 Calculation data of three principal stresses (MPa).

4.3 Natural vibration characteristic

Firstly, the model of drillstring combination is established. Apply boundary conditions to the model, drillstring has one end hinge branch, the other end movable hinge, the centralizer movable hinge. Take the previous order eigenvalues, the natural frequency of the BHA will be obtained as shown in Table 2.

Table 2 Natural frequency of the BHA.

As shown in Table 2, for a drillstring with a large length, the natural frequency or critical speed is relatively dense, and it is difficult to completely avoid the critical speed. At this time, it is necessary to avoid the critical speed as much as possible, for example, the intermediate value of the two-step critical speed can be taken. In the actual drilling process, the drilling tool speed should be far away from the critical speed. The assembly dangerous speed for natural frequency of 1, 2, 3 and 4 order are corresponding to the rotational speed, i.e. 30, 73, 128 and 196 r min⁻¹. In the work under the rotational speed of drilling tools, resonance is inclined to be caused, and then lead to drilling tool failure fracture. In drilling process, rotary speed needs to stay away from these dangerous values.

4.4 Dynamic analysis

When the speed is 60, 100 and 140 r min⁻¹ respectively, the changes of the maximum stress point in different cycles are shown in Fig. 5. Stress changing period and maximum stress differ with different rotating speed. When the speed is 140 r min⁻¹, the stress peak is the largest, up to 60 MPa. When the speed is 100 r min⁻¹, the stress peak is the minimum.

Figure 5 Dynamic stress time series of each node at different rotary speeds.

Additionally, when the rotating speed is 60, 100 and 140 r min⁻¹ respectively, the drill bit lateral force changes in different periods are shown in Fig. 6. When the speed is 140 r min⁻¹, the lateral force reaches a maximum of 40 kN at 6.95 s. At different speeds, the lateral forces tend to stabilize with time. When the rotating speed is 140 r min⁻¹, displacement of assembly at various points in a cycle is shown in Fig. 7a, due to the constraints on both ends of the assembly, there is no displacement in X and Z directions, the other drillstring point has larger displacement changes. When the speed is 140 r min⁻¹, as shown in Fig. 7b, the stress at each node of the drillstring constantly changes within a period. The node with maximum stress is in the middle of the drillstring as indicated by the arrow. The maximum stress reached about 78 MPa at 5.863 s.

Figure 6 Lateral force time series of each node.

Figure 7 The variation of the displacement and the stress of each node in the drill assembly with time.

4.5 Safety evaluation

When the rotary speed is 140 r min⁻¹, the dangerous stress point is 14.62 m away from the bottom of the drillstring. The change of stress over time is shown in Fig. 8.

As shown in the Fig. 8, the maximum stress ${\mathit{\sigma }}_{max}=\mathrm{60.8264}$ MPa, the minimum stress ${\mathit{\sigma }}_{min}=\mathrm{0.4207}$ MPa, and equivalent stress amplitude is ${\mathit{\sigma }}_{\mathrm{da}}=\frac{{\mathit{\sigma }}_{max}-{\mathit{\sigma }}_{min}}{\mathrm{2}}=\frac{\mathrm{60.8264}-\mathrm{0.4207}}{\mathrm{2}}=\mathrm{30.20}$ (MPa). Since the σ is obviously less than σ_s and σ_da<σ_a, it meets the fatigue strength requirement.

Figure 8 Curve of dynamic stress over time of each node at rotary speed is 140 r min⁻¹.

4.6 Experimental vibration data analysis

Direct measurement of stress values is difficult, besides the theoretical research, the measurements of drillstring vibration in in field were carried out by drilling engineer. To understand the intensity of drillstring vibration, Alley and Sutherland (1991) and Tian (2016) used Measurement While Drilling (MWD) to measure vibration parameters, the results show that the drillstring vibrates violently while drilling. The latest measurements were measured by the micro vibration logging tool, which was directly mounted on the drill bit to avoid interference of other downhole tools, so the vibration characteristics of drillstring and the fluctuation of WOB will be truly reflected by the measurements.

The rotary steerable system (RSS) with push-the-bit mode is installed near the bit for directional drilling, the whole drillstring of the drilling system is rotated from the surface by a top drive, which is typically a hydraulically driven mechanism. During directional drilling using RSS, specialized downhole equipment is employed to replace conventional directional tools such as mud motors. For RSS using push-the-bit mode, its push pads in implementing agency alternately in turn pushed against the borehole wall, making the bottom hole constantly under a cycle of nonlinear damping force, by which the movement of the RSS would be transformed to chaos and disorder. Thus, it is necessary and urgent for us to do some thorough research about the dynamic state of RSS. High-sample-rate downhole dynamics sensors were placed at 0.46 m above the bit in the BHA.

We develop a strap-down measurement while drilling (MWD) surveying system that incorporates three-axis magnetometers and three-axis accelerometers arranged in three mutually orthogonal directions (Russel and Russel, 1979). The sensors are installed inside nonmagnetic drill collar that can avoid the external magnetic interferences (Rehm and Garcia, 1989). a_x, a_y, a_z are defined as survey signals of triaxial accelerometers on the x, y, z axis respectively. m_x, m_y, m_z are defined as survey signals of triaxial magnetometers on the x, y, z axis respectively, the sampling frequency f_s is 1 kHz. Assume that the Earth's magnetic field strength as M. Obviously, $M=\sqrt{m_x^{\mathrm{2}}+m_y^{\mathrm{2}}+m_z^{\mathrm{2}}}$. Under certain sample frequency, measuring signal is time series and can be expressed as time function. We can use $a_h=\pm \sqrt{a_x^{\mathrm{2}}+a_y^{\mathrm{2}}}$ define the lateral vibration of the BHA and use a_z express the longitudinal vibration.

Figure 9 The time series of drillstring vibration and the spectrogram of x and y signals of accelerometers, in the case of safety drilling.

Figure 10 The actual photo of the drillstring damage (a) and the vibration signal measured in the early stage of the drillstring damage (b).

As shown in the Fig. 9, the time series of lateral vibration and longitudinal vibration of drillstring were illustrated, and the spectrogram of x and y signals of accelerometers were illustrated subsequently. The vibration data is measured in the case of safety drilling. In the secondary development software as mentioned above, the σ is less than σ_s and σ_da<σ_a, indicates that the working condition is safe. In the other case, if the vibration state of the drillstring is unsafe, the measured vibration signals level will increase, as shown in the Fig. 10, the vibration signal gradually increases until the drillstring is damaged. Simultaneously, we can calculate by the software we developed in this paper, σ≥σ_s or ${\mathit{\sigma }}_{\mathrm{da}}\le {\mathit{\sigma }}_{-\mathrm{1}}$, reveal the vibration intensity of drillstring is too high to safe.

5 Conclusions

The vibrations induced in the drillstring are the main reasons for the drillstring failure. To develop practical and commercial computing software about drillstring dynamics is difficult, the model established must take lithology, bit type, drill assembly, damping and contact into consideration, meanwhile through a lot of experiments and experiment data analysis. In this paper, we developed a secondary development software based on commercial finite element software, which reduces large amount of work for calculation method and program verification. The dynamic characteristics of the inherent frequency and dynamic stress of the bottom hole assembly (BHA) were calculated and analyzed, and risk assessment method based on the quantitative vibration intensity was established. Considering the contact and collision between drillstring and borehole wall, taking the whole drillstring assembly system as the research object, the dynamic equation is solved by explicit central difference method in Abaqus, and the dynamic response of drillstring system is obtained. Then a risk assessment method based on quantitative vibration intensity was established and dynamic characteristics and safety of typical drillstring were analyzed. Using the secondary development function of Abaqus to develop specific graphical user interface and command flow to simplify, it is simpler and more effective to analyze specific engineering problems.

The “danger assessment method based on quantitative vibration intensity” established in this paper firstly requires avoiding the critical speed and selecting the ideal recommended speed; if it is not fully realized, the dynamic stress is evaluated, and the vibration intensity is determined according to the fatigue allowable stress and the yield limit of the drillstring. Divided into three grades to assess the reliability of the drill assembly with dynamic stress levels. This method combines the point control to avoid the resonance speed and the interval control based on the dynamic stress level, and forms the quantitative evaluation principle of the vibration risk of the drillstring.

The vibration and internal stress of the drillstring have a positive correlation, downhole measurement can only get the vibration signals and cannot get the magnitude of the stress value, simultaneously, internal stress is easily obtained by theoretical calculation. In this paper we developed a software for secondary development that can quickly get stress values. The combination of the two methods, ability to provide an optimal solution for the safety assessment of drilling tools. In addition, the effectiveness of our reliability assessment method is also proved by experimental vibration data.