Dynamic Model and Quantitative Analysis of Stick-Slip Vibration in Horizontal Well

Stick-slip is very harmful to the service life of drillstring. ,e extended Hamilton principle is applied in the paper. ,en, finite element method (FEM) is employed to describe the model. ,e drillstring-borehole impact and friction, fluid-structure interaction, bit-rock interaction, and gravity are considered in this model. ,e influence of axial and torsional excitation on stickslip is analyzed. ,e nonlinear motion predicted by the model is consistent with the observation results in the experiments. ,e research shows that the fluctuation amplitude of the bit angular velocity also increases along with the increase of driving angular velocity (torsional excitation). However, both the ratio of the maximum angular velocity of the stick-slip vibration and the fluctuation of the angular velocity are continuously reduced. Meanwhile, the strength of the stick-slip vibration has a tendency to slow down. As the axial load (axial excitation) increases, the fluctuation of the maximum angular speed of the stick-slip vibration does not change significantly, but the smaller load causes a smaller stick duration.


Introduction
Rotary drilling system for exploration and development of oil and gas reservoirs in petroleum engineering is widely applied. e drillstring is a series of drilling tools, which may include drill pipes, drill collars, stabilizers, and auxiliary tools. e axial load, torque, bit-rock interaction, slurrydrillstring interaction, and slurry damping have negative effects on the drillstring. Due to the complex environment, a statistic data will produce some unwanted vibrations, in which stick-slip is one of the most serious ones [1][2][3][4][5]. e drillstring often experiences severe torsional vibration and stick-slip vibration, in which the angular velocity near the bit can be increased to two or more times of the applied angular velocity at the top of the drillstring [6][7][8].
Belokobyl'skii and Prokopov research the friction that induces self-excited stick-slip in the early days [9]. ey present a single degree of freedom system with dry friction and think that the nonlinear friction can cause stick-slip. Rappold also analyzes a self-excited oscillating system under dry friction and studies the stick-slip vibration with chaos [10]. A nonregularized dry friction model of bit-rock interaction is applied to analyze the stick-slip motion by Hugo. e mean deviation of angular velocity is obtained. en the stability regions are identified [11]. Navarro-López and Licéaga-Castro develop a kind of control method on dynamical sliding mode to suppress stick-slip vibration in vertical well. And the mathematical model with four DOFs is presented. Furthermore, they present a new procedure to analyze stick-slip vibration based on stability characteristics and the relationship of various sliding motions [12]. Tang et al. investigate the influence of static and kinetic friction on stick-slip vibration. And a lumped-parameter model is developed [13]. Sarker et al. develop a bond graph dynamic model in horizontal well. e friction in this model is divided into two parts, frictions in deflecting section and horizontal section. Coupled behavior among bit-rock interaction, angular velocity of bit, and weight on bit (WOT) can be predicted in this model [14].
In some researches, slurry is taken into account in the models [15][16][17]. Ertas et al. develop a general dynamical model of drillstring taking into account the viscosity of slurry [15]. Yigit and Christoforou develop a mathematical model to obtain the dynamic behavior with the influence of drilling slurry. e results show that an additional axial excitation can effectively suppress stick-slip vibration [16]. Ytrehus et al. propose that the friction between drillstring and borehole wall is dependent on the physical property of drilling slurry [17].
Bit-rock interaction is also an important part in the horizontal well [4,[18][19][20]. Kovalyshen applies drag-bit model to study coupled axial-torsional vibration. e interaction between rock and bit includes cutting and friction. Axial and torsional vibrations are coupled through bit-rock interaction [21]. Kamel and Yigit propose a more realistic nonlinear model to describe cutting and contact components with axial and torsional motions [1]. A drag-bit model with two degrees of freedom is presented by Zamanian et al. in order to study the stick-slip vibration of drillstring. e research shows that bit-rock interaction has an effect on stick-slip vibration, which can be avoided by optimizing bit parameters [22]. Kovalyshen thinks that the stick-slip vibration is caused by cutting motion and wear between rock and bit.
In the present researches, most of them are aimed at double coupled vibration or vertical well. However, there are few studies on fully coupled vibration of horizontal well with fluid, friction, and bit-rock interaction. In this paper, the nonlinear model on fully coupled vibration in horizontal well under various complex boundary conditions is applied. And the effects of axial excitation and torsional excitation on stick-slip vibration are researched. Finally, the mathematical model is verified by experiments.

Mathematical Model
Because the drillstring in horizontal well is subjected to shear deformation, Timoshenko beam theory is employed in this paper. Moreover, the drillstring is also subjected to the friction from borehole wall, drilling slurry, loading at the end, and the driving torque at the top, which make the mathematical model of drillstring more complicated. erefore, the motion equation is deduced by an energy method-extended Hamilton principle, which is of great advantages in dynamic modeling with complex boundary conditions [23].

Basic Model.
e model presented in this paper includes three kinds of vibrations including axial vibration, torsional vibration, and lateral vibration. ey are shown in Figure 1.
According to the characteristics and the purpose of model, it is a crucial step in modeling to make necessary and reasonable simplifications and assumptions. If all the factors are taken into account, the whole model will be too complicated. us, the consumption of time and resources are huge. And the improvement of accuracy is not so much because of the complexity of the model. e following assumptions are made in this paper with a certain calculation accuracy: (1) e drillstring is horizontal and the cross section is a regular circle (2) e material of drillstring is isotropic and in linear elasticity (3) e drilling slurry is ideal

Motion Equation.
e extended Hamilton principle is used to derive the motion equation of drillstring. e extended Hamilton principle expression is where W is the work done by the external forces, T is the kinetic energy, and U is the strain energy. e model is discretized by finite element method. e drillstring has a total length of L and is equally divided into n sections that has a length of l. ere are two nodes and 12 degrees of freedom in each unit that is shown in Figure 2.
In Figure 2, u, v, w, and θ denote the displacement of axial, torsional, lateral, and angular displacement, respectively.
Make its local coordinates be ξ � x/l e . Because the model is discretized, the element displacement can be given as where t denotes time. N u , N v , N w , N θ x , N θ y , and N θ z denote the shape functions of each displacement, respectively.

Work Done by Gravity.
In the process of modeling the drillstring, gravity is a conservative force that largely affects its dynamic behavior. Work done by gravity is not related to the path of the object, but only to its initial and final positions. e work done by the conservative force gravity can be written as where A is the cross-sectional area of drillstring, g denotes acceleration of gravity, and ρ denotes the drillstring density. Finite element method and variation are used to discrete each element. e force vector can be written as where N T v is the transposed matrix of the shape function N v .

Kinetic Energy of
Drillstring. e kinetic energy expression of drillstring can be given as

Shock and Vibration
where I is cross-sectional moment matrix of inertia, r is angular velocity vector, and s is velocity vector. e expression can be written as I � where I 1 and I 2 are polar moment of inertia and crosssectional moment of inertia, respectively. Finally, the kinetic energy of drillstring in variational form can be written as where δu, δv, δw, δθ x , δθ y , and δθ z denote the variation of u, v, w θ x , θ y , and θ z , respectively. e second derivative term is the mass matrix of the model and the rest is the force vector:

Strain Energy of Drillstring.
e drillstring is a system of mass points whose positions have changed after deformation.
us, the strain energy of the whole system has changed. According to the principle of elastic mechanics, the strain energy can be written as where σ denotes stress tensor, ε denotes strain tensor, and V denotes drillstring volume.

Friction and Impact.
When the drillstring rotates in a horizontal well, it will bend in some cases. us, friction and impact may occur between the drillstring and the borehole. In this paper, elastic force is employed to describe the friction between drillstring and borehole wall [24]. e elastic forces in the y and z directions can be given as Shock and Vibration 3 where k im is impact factor, Ψ im (r) is coefficient of determination, R e denotes external radius of drillstring, R h denotes borehole radius, and r is the radial displacement of drillstring. As shown in Figure 3. e friction between the borehole and the drillstring during rotation is modeled by friction torque in x direction: where μ is the coefficient of kinetic friction, sign denotes symbolic function, F n denotes normal contact force, and e variation form of the work done by the impact can be written as where x n denotes coordinate at node n and q e denotes displacement vector.

Bit-Rock Interaction.
Bit-rock interaction is composed of the torque on bit and the weight on bit. ey can be decomposed into the forces transmitted by the cutting interface (the subscript c) and the wear flats (the subscript f ) [25]. us, they can be written as T c and W c can be given as , B n denotes the number of blade of bit, σ denotes maximum contact force, ζ denotes the cutting force orientation description, ε denotes intrinsic energy, and a denotes the bit radius. e interaction can be composed of four parts: is section dictates that the minimum forces (W and T) are applied: us, the variational form of the work done by cutting and wear is In Figure 3, F N and F A denote viscous hydrodynamic force and lateral nonviscous hydrodynamic force, respectively. F it and F in denote the tangential and normal component: where M f is the mass per unit, A i � (D 2 i /4)π, and p i is the pressure of internal flow. e force caused by flow outside the drillstring along x direction can be written as where

Shock and Vibration
e resultant force from internal flow in z direction can be written as e resultant force caused by external flow along z direction is where , ρ f denotes drilling slurry density, f z denotes frictional viscous force, and U ex and U i denote the flow velocity of the external and internal slurry. FEM is employed to discretize the physical model. en, the second derivative term of displacement is considered as mass matrix: e first derivative term is considered as the damping matrix: e first term is considered as the element stiffness matrix: e rest of the term is the force vector:

Discretization of Motion Equation.
e motion equation can be written as where M, C, K, and f denote the mass matrix, damping matrix, stiffness matrix, and force vector, respectively. Each matrix can be written as follows: A damping model is presented in this paper to estimate the dissipation of system. A comprehensive estimation of the energy dissipative matrix can be written as [26] where Θ 1 and Θ 2 are the additional damping coefficient.

Numerical Results and Discussion
ere may be many factors that affect stick-slip vibration, such as the angular velocity, axial force, friction, and the load of bit.
e numerical simulation about the influence of angular velocity and axial force on stick-slip vibration of the drillstring is conducted in this paper. Figure 4: External flow forces.

Relationship between Torsional Excitation and Stick-Slip
Vibration. Because FEM is used for discretization, the six degrees of freedom at the last node of drillstring can be obtained. en, the angular velocity is extracted to obtain the angular velocity curve in time domain, just as shown in Figure 6.
As shown in Figure 6, the simulation speeds are 17.9 r/ min, 35.7 r/min, 53.6 r/min, 74.4 r/min, 89.3 r/min, and 107.1 r/min, respectively. When the driving speed is at the low-speed range (17.9 r/min; 35.7 r/min), the regularity of stick-slip vibration is poor. When the speed is 17.9 r/min, the angular velocity of the bit fluctuates greatly in the initial stage, which causes a large dynamic load and the premature failure of drillstring and bit.
When the velocity is 35.7 r/min, the regularity of the bit is the worst and the angular velocity fluctuates violently. When the driving velocity is in the middle section (53.6 r/ min; 74.4 r/min), the angular velocity of the bit fluctuates softly, which is more regular than the driving velocity is at low speed. When the driving speed is at the high-speed range (89.3 r/min, 107.1 r/min), the amplitude of the angular velocity drops off. In particular, when the driving velocity is 89.3 r/min, the fluctuation of the bit angular velocity is small. us, it is an ideal driving velocity. In order to observe the relationship between the maximum angular velocity fluctuation of the bit and the angular velocity caused by stick-slip vibration more clearly, the maximum amplitude of numerical simulation under different driving velocities is obtained and compared with the different driving velocity, just as shown in Figure 7.
Driving velocity and the maximal angular velocity of the numerical simulation are shown in Figure 6. With the drive speed increases, the maximum amplitude of the angular velocity also increases. However, the proportion between the maximum angular velocity and the driving velocity is decreasing, whose ratios are 4.8, 4.0, 3.7, 2.6, 2.43, and 2.4, respectively. It can be seen that the angular velocity fluctuation of the bit caused by stick-slip vibration constantly decreases, and the intensity of stick-slip vibration tends to be weakened. e speed fluctuation caused by stick-slip vibration tends to be modest at the high-speed range. us, a higher driving angular velocity can suppress the strength of stick-slip vibration.

Relationship between Axial Excitation and Stick-Slip
Vibration. Axial excitation mainly refers to the axial force exerted on the top of the drillstring. In this section, the axial excitation will be exerted by changing the axial load on the top of the drillstring to study the influence of different axial excitation on the dynamic behavior of the drillstring.
Since there is no independent load term in the motion equation, an initial equivalent axial force is added when calculating the force vector of element. Finally, the angular velocity of the bit is obtained.
As can be seen from Figure 8, the fluctuation of the maximum angular velocity of the stick-slip vibration does not change significantly with the increase of the axial load. erefore, it can be seen that the axial force does not affect the maximum of angular velocity fluctuation and does not have too much influence on the drillstring. However, in the global stick-slip vibration, the frequency of stick-slip vibration caused by smaller axial load is lower than that of higher axial load. Moreover, the stick state caused by the bigger axial load lasts longer than the smaller load.

Experimental Setup.
Similarity theory is applied to design an experimental apparatus in order to verify the mathematical model [27]. Schematic diagram of experiment apparatus of drillstring constrained in a horizontal pipe is shown in Figure 8. Based on the similarity theory and the functions of the test equipment, the drillstring dynamic behavior test system mainly includes five parts: driving system, pipe system, testing system, cycling system, and loading system, just as shown in Figures 9 and 10. e similarity criteria are deduced by the dimension analysis method. e dimensions of the parameters applied in this paper are shown in Table 1. And the dimension matrix of the system is shown in Table 2.
e dimension matrix of the system can also be written as some equations just as shown in a 2 + a 3 + a 4 + a 6 � 0, e following equations of π are obtained on the basis of π matrix: According to the similarity theory, there are two kinds of similarity ratios. One is the geometric similarity ratio. e other one is the physical similarity ratio. All physical quantities with unit of length in nonambiguity take the same geometric ratio, such as the inside and outside diameter and length of the drillstring. us, the length similarity ratio c_l � l/l � 1: 10, where l and l denote the length in simulation and experiment, respectively. e parameter with an angular unit in nonambiguity remains the same. In this paper, the similarity ratio of angle c_α � α/α � 1, where α and α denote the angle used in simulation and experiment, respectively.
Based on the geometric similarity ratio, the parameters used in simulation are given as follows: In order to ensure that the test device is consistent with the real working condition of the drillstring, it is necessary to determine the main body of the experiment apparatus; the material of the drillstring in experiment and real material in simulation (real working condition) have the same or similar stress-strain curves in order to make sure of the accuracy of the results. us, the models of simulation and experiment   Shock and Vibration are similar in physics in this paper. e density similarity ratio c_ρ � ρ/ρ � 0.14, where ρ and ρ denote the density of drillstring in simulation and experiment, respectively. e similarity ration of elasticity modulus c_E � E/E � 0.011, where E and E are the elasticity modulus in simulation and experiment, respectively. Based on the physical similarity ratio, the parameters used in simulation are given as follows: Shock and Vibration e density of drillstring, ρ � 7.85 g/cm 3 e elasticity modulus of drillstring, E � 210 GPa e parameters used in experiment are given as follows: e density of drillstring, ρ � 1.10 g/cm 3 e elasticity modulus of drillstring, E � 2.30 GPa e other similarity ratios in this paper are c_ω � ω/ω and c_P � P/P. P and ω are the axial load and angular speed in the experiment, respectively. P and ω are the axial load and angular speed in the numerical simulation, respectively. In the paper, c_ω � 2.8 and c_P � 1/9130. e experimental and the simulation parameters, applied in this paper, are shown in Tables 3 and 4, which satisfied the similarity ratios. e function of the loading system (Figure 11(a)) is to simulate the contact between the bit and the bottom hole. And a pressure sensor is installed at the bottom to measure the axial force acting on the bit in real time. e driving system is shown in Figure 11(b). e servo motor increases the torque and reduces the rotating speed through the planetary reducer and then transmits the rotating power through the torque sensor to drive the drillstring, which can simulate the rotation exerted on the drillstring in the drilling process of a horizontal well. In the experiment, the rotating speed of the servo motor can be controlled by adjusting the servo driver; different rotating speeds are applied in the horizontal drillstring test device. e influence of the rotating speed on the dynamic state of the drillstring can be detected. In addition, the equipment for applying WOB is also integrated into the system to achieve the purpose of convenient operation.
In order to detect the dynamic state of drillstring, 5 detection points are set in the experiment apparatus. ey make up the testing system (Figure 11(c)). Due to the presence of water in the annular space, it is necessary to select eddy current displacement sensors that can be waterproof and realize noncontact measurement. e test points include a total of 5 (15 eddy current displacement sensors): 3 sensors (in x, y, and z directions, respectively) are installed orthogonally on each test point to, respectively, detect the displacement of different cross sections of the drillstring in different directions.
In the process of drilling, the existence of circulating slurry will also have a certain impact on the dynamic behavior of the drillstring.
e circulatory system (Figure 11(d)) can simulate friction slurry at different flow speed and the effect of slurry on the dynamic state of drillstring. Water is used to simulate slurry in this paper. e circulatory system includes a tank, a pump, and a flowmeter. e drillstring and the plexiglass tube form the channel of slurry. e increase of flow and flow rate of slurry in the drillstring can be observed in the flowmeter through adjusting the pump. e workflow of the system is as follows: pump the slurry into the drillstring through the flow meter, slurry flows into the annular space between the drillstring and the plexiglass tube, slurry returns to the tank through ball valve, and complete the slurry circulation.

Relationship between Torsional Excitation and Stick-Slip
Vibration. In order to verify the model, an experiment is carried out about the influence of the driving velocity on the stick-slip vibration according to the simulation results of the mathematical model. e angular velocity of the bit is measured by a sensor at the bottom of the drillstring. e curve is shown in Figure 12.
Experimental angular velocity is 100 r/min, 200 r/min, and 300 r/min, respectively. It can be observed that when the driving velocity is at the low-speed range, the angular velocity of the bit fluctuates violently. However, when the driving velocity is 200 r/min, although the amplitude fluctuation is large, the ratio of the maximum fluctuation to the driving speed is still small. When the driving speed is 300 r/ min, the test speed rises and fluctuates regularly. e driving speed, the fluctuation of the test speed, and the relationship between them are shown in Figure 13 in order to identify the influence of the driving velocity on stick-slip vibration. e comparison between the experiment and the simulation is shown in Figure 13. e numerical results in the brackets are the driving velocity of the numerical simulation under the similarity theory. It can be found that, with the driving velocity increases, the maximum angular velocity in experiment and simulation will also increase. A larger driving velocity at the top of the drillstring will inevitably lead to a larger rotating speed of the bit. However, the ratio between them is decreasing. It can be considered that the angular velocity fluctuation caused by stick-slip is decreasing and closer to the driving angular velocity. And the detriment to the drillstring and bit is reduced. e numerical simulation is consistent with the experimental results. us, the result proves the correctness of the mathematical model.

Relationship between Axial Excitation and Stick-Slip
Vibration. According to the structure of the experiment apparatus, the dynamic behavior of drillstring under different axial excitations is studied through changing the axial force, just as shown in Figure 10. Figure 14 shows the angular velocity of the bit with different loads. e angular velocity curve is almost matched when the axial load is 0.5 kgf∼2.5 kgf. And there are no significant changes on the maximum of the angular velocity (only short-term fluctuation occurs at the initial stage). e stick-slip strength of the bit cannot be clearly distinguished under the three conditions. erefore, the total stick-slip time is applied to measure the intensity of the stick-slip. Figure 15 shows the stick-slip duration under simulation and experiment.
According to statistic data, the total stick durations of the bit are 11.5 s, 13.6 s, and 14.3 s under different axial load, respectively. It can be found that, with the initial axial load     increase, the total stick duration shows an upward trend. It also indicates that the strength of stick-slip is increasing. e numerical simulation results are performed at different driving velocities according to the similarity theory. e trends of simulation and experiment are exactly the same. us, the experimental results are consistent with the numerical simulation. Experiments have proved the correctness of the numerical simulation.

Conclusions
e dynamic behavior of the drillstring under axial and torsional excitation is researched in this paper. e mathematical model is discretized by FEM. A numerical simulation has been conducted in MATLAB and verified by the experiment. e following conclusions can be obtained: (1) With the drive velocity increases, the maximum of the numerical simulation angular velocity also increases, but the ratio of the maximum stick-slip angular velocity to the driving angular velocity decreases. e intensity of stick-slip vibration tends to be weak. And the velocity fluctuation is modest at the high speed. (2) e axial load does not affect the maximum of angular velocity fluctuation. However, in the global stick-slip vibration, the stick-slip frequency of drillstring caused by smaller axial load is lower than that of higher axial load. Moreover, the stick state caused by the larger axial load lasts longer than that at the low-speed range. (3) e design of the experimental apparatus in this paper is feasible. And the experimental results are in good agreement with the simulation results. Finally, the mathematical model is proved to be correct.

Data Availability
Data are not made available as the NSF supporting the project is confidential.

Conflicts of Interest
e authors declare that they have no conflicts of interest.