Design of Telerobotic Drilling Control System with Haptic Feedback

Abstract

The paper deals with the design of control algorithms for virtual reality based telerobotic system with haptic feedback that allows for the remote control of the vertical drilling operation. The human operator controls the vertical penetration velocity using a haptic device while simultaneously receiving the haptic feedback from the locally implemented virtual environment. The virtual environment is rendered as a virtual spring with stiffness updated based on the estimate of the stiffness of the rock currently being cut. Based on the existing mathematical models of drill string/drive systems and rock cutting/penetration process, a robust servo controller is designed which guarantees the tracking of the reference vertical penetration velocity of the drill bit. A scheme for on-line estimation of the rock intrinsic specific energy is implemented. Simulations of the proposed control and parameter estimation algorithms have been conducted; consequently, the overall telerobotic drilling system with a human operator controlling the process using PHANTOM Omni haptic device is tested experimentally, where the drilling process is simulated in real time in virtual environment.

1. Introduction

Drilling a borehole is a common method for extracting oil, gas, and natural resources from beneath the surface of the earth. Conventional oil well drilling has made significant progress over recent years, and currently is one of the most automated processes in the oil and gas industry. However, there are still some fundamental challenges associated with the drilling. One of the challenges is the choice of vertical penetration velocity of the drill bit. For efficient drilling operation, this velocity must depend upon the type of rock beds drilled. In particular, the velocity must be adjusted when mechanical characteristics of rock strata change. Often, it is difficult to estimate in real time the relative position of the drill bit with respect to different rock layers and, therefore, hard to predict the mechanical characteristics of the rock formations.

The goal of this research is to design a telerobotic system with haptic feedback for control of the drilling process. Telerobotics for drilling well is a relatively novel idea, and it is substantial endeavor to automate one of the fundamental processes in the extraction of energy and resources. As telerobotics is integrated with drilling, it can greatly decrease the number of people working and monitoring operation on the site. This, in particular, can reduce the work site hazards. Also, telerobotics can bring actual analysis of in situ conditions (underground drilling environment) in real time to the human operator that works remotely, where (s)he will be able to monitor the current drilling conditions and, in particular, promptly enforce changes in the vertical speed of penetration of the drill bit in the oil well. Real-time control and optimization of the drilling speed are crucial for today's drilling industry, as it can reduce time and immense cost associated with the drilling an oil well. Introduction of haptic feedback would allow the human operator to feel the changes in mechanical characteristics of the rock and adjust the vertical velocity of penetration accordingly.

In this paper, we address the problem of design of control algorithms for virtual reality based telerobotic system with haptic feedback that allows for the remote control of the vertical drilling operation. Based on a simplified mathematical model of the drilling process, control algorithms are designed which allow to achieve a desired rate of the vertical penetration, regardless of the mechanical properties of the rock. The control design includes an online parameter estimator of the intrinsic specific energy which is a parameter that describes the hardness of the rock. All these algorithms are consequently used in the design of a telerobotic drilling system with virtual environment-based haptic feedback that allows the human operator to feel the stiffness of the rock in contact with the drill bit. Simulations and semiexperimental results are performed which confirm the validity of the theoretical developments.

The potential application domain of this research is not limited to onshore/offshore oil well drilling, but the same principles can be applied, in particular, to different types of mining robots [1], telerobotic systems for dredging and mining ocean [2–5], surgical drilling [6], and telerobotic systems for drilling the extraterrestrial terrain to discover and research the minerals and composition beneath [5, 7].

The structure of the paper is as follows. In Section 2, a mathematical model of the drilling process is derived which is subsequently used for the control design. Section 3 deals with the design of control algorithms for rotational and translational motion of the drilling systems, as well as the design of an online parameter estimator of the intrinsic specific energy of the rock. In Section 4, the structure of a telerobotic drilling system is described and the corresponding experimental results are presented. Finally, in Section 5, some conclusions are given and possible future directions are formulated.

2. Mathematical Model of Drilling System

In this section, mathematical models that describe the drilling system are presented. Specifically, the mathematical model of drill string and drive system is described in Section 2.1, while the model of rock cutting and penetration is the subject of Section 2.2.

2.1. Mathematical Model of the Drill String and Drive System

The drill string is the assembly of rotating pipes which are responsible for transmitting rotation and weight to the bit and bridge up a connection between the bottom hole tools [8]. The components of a drill string along with drill pipes and the bottom hole assembly (BHA) are shown in Figure 1. A number of simplified mathematical models for drill string and drive systems were proposed in the literature, such as [9–12]. The model used in our work was developed in [9]. This model describes the drill string as a simple torsional pendulum, where the drill pipes are represented as torsional springs and the bottom hole assembly is described as a rigid body with inertia. The model is based on the following simplifying assumptions.(1)The bottom hole assembly and the drill bits behave like rigid bodies.(2)The moment of inertia of the drill pipe is considered to be small in comparison with the moments of inertia of the bottom hole assembly and the rotary table and, therefore, neglected.(3)The nonzero time propagation of the torsional force disturbances along the drill string is neglected. The forces assume to propagate instantaneously along the drill string.

Figure 1 

Drill string components [8].

Under the above described assumptions 1–3, the whole drill string and drive system with equivalent electro-mechanical components can be represented by its structural diagram shown in Figure 2. This system is described by the following mathematical model [9]. First, the motion of the drill string is described by the following equation: 𝐽₁𝜙̈₁+𝑐₁𝜙̇₁+𝑘𝜙₁−𝜙₂−𝑇=0. (1) Here, 𝜙₁ is the angular displacement of bit and drill collars (BHA), 𝜙₂ is the angular displacement of the rotary table, 𝐽₁ is the equivalent moment of inertia of the collars (BHA) and the drill pipes, coefficient 𝑐₁ represents equivalent viscous damping, 𝑘 is the equivalent torsional stiffness of the drill pipes, and 𝑇 is the torque-on-bit (TOB) generated during the rock cutting process (see Section 2.2 below). The dynamics of the rotary table and drive system is described by the following equation: 𝐽₂𝜙̈₂+𝑐₂𝜙̇₂−𝑘𝜙₁−𝜙₂−𝑛𝑇_𝑚=0, (2) where 𝐽₂ is combined moment of inertia of the rotary table and of the rotor of the electric motor coupled together with a gearbox that has 1∶𝑛 gear ratio, 𝑐₂ is aggregated damping of all the components of the drive system, and 𝑇_𝑚 is the motor torque. Finally, the electric motor is described by the following equations: 𝐿𝐼̇+𝑅𝐼+𝑉_𝑏−𝑉=0,𝑉_𝑏=𝐾𝜙̇₃=𝐾𝑛𝜙̇₂,𝑇_𝑚=𝐾𝐼, (3) where 𝐼 is the armature current, 𝐿 is an equivalent armature inductance, 𝑅 is an equivalent armature resistance, 𝑉_𝑏 is the back emf, 𝑉 is the armature voltage, 𝜙̇₃ is the rotor angular velocity, and 𝐾 is a constant that depends upon the motor characteristics.

Figure 2 

Representation of drill string/drive system with mechanical and electrical components [9].

By combining all the above equations, the complete drill string/drive system can be written in the following state space form: 𝜙̇₁𝜔̇₁𝜙̇₂𝜔̇₂𝐼̇=01000−𝑘𝐽₁−𝑐₁𝐽₁𝑘𝐽₁𝑘0000010𝐽₂0−𝑘𝐽₂−𝑐₂𝐽₂𝐾𝑛𝐽₂000−𝐾𝑛𝐿−𝑅𝐿𝜙₁𝜔₁𝜙₂𝜔₂𝐼+0−𝑇𝐽₁00𝑉𝐿. (4) Here, 𝜔₁ and 𝜔₂ are the angular velocities of the drill bit and the rotary table, respectively. Equation (4) is valid when the drill bit rotational velocity is greater than zero, that is, 𝜔₁>0. In order to reduce the number of equations, a variable 𝜙 is introduced as the difference of 𝜙₂ and 𝜙₁. In this case, the original system can be rewritten in the following reduced state space form: 𝜔̇₁𝜙̇𝜔̇₂𝐼̇=−𝑐₁𝐽₁𝑘𝐽₁000−1010−𝑘𝐽₂−𝑐₂𝐽₂𝐾𝑛𝐽₂00−𝐾𝑛𝐿−𝑅𝐿𝜔₁𝜙𝜔₂𝐼+−𝑇𝐽₁00𝑉𝐿. (5) Equation (5) defines the reduced order model of the drill string and drive system. The model (5) is used for the control design below.

2.2. Rock Cutting and Vertical Penetration Models

A standard drill bit usually exhibits two kinds of motions: rotational along its axis of rotation and vertical motion while penetrating through the rocks. As described in [13], in the normal mode of operation of the drill bit, the bit rotational velocity 𝜔 is parallel to its axis of rotation, and the penetration velocity 𝑣 is directed vertically straight through the rocks. Similarly, the weight-on-bit 𝑊 acts in the vertical direction and the torque-on-bit 𝑇 is applied in parallel to the direction of rotation of drill bit. The cutting components of the weight-on-bit and torque-on-bit depend on the radius of PDC drill bit 𝑎, intrinsic specific energy 𝜖, a parameter 𝜁>0 which represents the ratio of the vertical force to the horizontal force between rock and cutter contact surfaces, and the depth of cut 𝑑. The depth of cut 𝑑 plays significant role in the equations to follow that describe the cutting components of the torque-on-bit 𝑇 and the weight-on-bit 𝑊. The equations for these two cutting components are as follows [13]: 𝑇^𝑐=12𝑎²𝑊𝜖𝑑,_𝑐=𝑎𝜁𝜖𝑑. (6)(7)

In this work, the system is developed under simplifying assumption that the friction effects are negligible. In this case, both variables 𝑇≈𝑇^𝑐 and 𝑊≈𝑊^𝑐 are proportional to the depth of cut 𝑑, according to (6) and (7). As illustrated in Figure 3, the depth of cut 𝑑 is the thickness of rock ridge in front of the blade. It is assumed that the drill bit has 𝑛 number of identical blades, and the difference of angular positions of these two successive blades is (2𝜋/𝑛). In this case, 𝑑 is the combined depth of cut of all 𝑛 blades in each revolution of drill bit, according to the formula 𝑑(𝑡)∶=𝑛𝑑_𝑛(𝑡), (8) where 𝑑_𝑛 is the depth of cut of each blade. The depth of cut for each blade is in turn defined according to the formula 𝑑_𝑛(𝑡)∶=𝑈(𝑡)−𝑈𝑡−𝑡_𝑛, (9) where 𝑈(𝑡) and 𝑈(𝑡−𝑡_𝑛) are the vertical positions of the drill bit at current time instant 𝑡 and a certain previous instant 𝑡−𝑡_𝑛, respectively [10, 11]. The delay 𝑡_𝑛 in the above formula is exactly the time that is required for the drill bit to rotate by an angle 2𝜋/𝑛 to achieve its current angular position 𝜙₁(𝑡); in other words, it also satisfies the following equation: 𝜙(𝑡)−𝜙𝑡−𝑡_𝑛=2𝜋𝑛. (10)

Figure 3 

Section of the bottom hole profile located between two successive blades [10].

Using (9) and (10) for calculating 𝑑(𝑡) would significantly complicate the control design. In this work, we simplify this problem by assuming that both the vertical and angular velocities change slowly; specifically, it is assumed that both 𝑣(𝜏)≡𝑈̇(𝜏) and 𝜔₁(𝜏)≡𝜙̇₁(𝜏) are approximately constant during each period 𝜏∈[𝑡−𝑡_𝑛,𝑡]. Using this assumptions, (9) and (10) can be rewritten as follows: 𝑑(𝑡)≈𝑛⋅𝑣(𝑡)⋅𝑡_𝑛,𝜔₁(𝑡)⋅𝑡_𝑛≈2𝜋𝑛. (11)(12) Combining (11), (12), and assuming 𝜔₁(𝑡)≠0, one gets the following approximate expression for 𝑑(𝑡): 𝑑(𝑡)≈2𝜋⋅𝑡(𝑡)𝜔₁(𝑡). (13)

The above formula has a singularity at 𝜔₁(𝑡)=0. To remove this singularity, note that the drilling occurs when both 𝜔₁(𝑡)>0 and 𝑣(𝑡)>0. On the contrary, 𝜔₁(𝑡)≤0, the drill bits do not cut the rock and therefore 𝑑(𝑡)≡0 in this case. Based on the above considerations, one can approximately define the depth of cut according to the formula 𝑑(𝑡)≈2𝜋⋅𝑣(𝑡)max𝜔₁(𝑡),𝜖₀, (14) where 𝜖₀>0 is sufficiently small positive constant. The formula (14) does not have singularity at 𝜔₁(𝑡)=0; it will be occasionally used for calculations of 𝑑(𝑡) instead of (13) in the cases where avoiding singularity is important (in simulations, etc.).

Finally, the vertical motion of the drill bit is described by the following equation [12]: 𝑀𝑑𝑣𝑑𝑡=𝑊_𝑠−𝑊−𝐻₀−𝐾_𝑓𝑣. (15) Here, 𝑣 is the vertical velocity of the drill bit, 𝑀 is the combined mass of the drill string and BHA, 𝐻₀ is the constant upward force applied from the top of drilling rig, and 𝑊_𝑠 is the submerged weight of the drill string and Bottom Hole Assembly (BHA). In this model, it is assumed that 𝑊_𝑠 and 𝐻₀ to be constants and defined their difference with another constant 𝑊₀ such that 𝑊₀=𝑊_𝑠−𝐻₀. Also, 𝑊 is the applied weight on bit from the interaction of rock defined by (7), and 𝐾_𝑓>0 is the coefficient of viscous friction.

3. Controller Design

The block diagram of the overall drilling system is shown in Figure 4. As it can be seen from this figure, the block diagram has a complex structure and consists of several interconnected subsystems. Specifically, the vertical motion subsystem is described by (15); the output of this subsystem is the vertical velocity of penetration 𝑣(𝑡). The subsystem that represents the rotational motion is described by (5); this subsystem has one control input which is the armature voltage 𝑉(𝑡) and one output which is the angular velocity of the drill bits 𝜔₁(𝑡). Both 𝑣(𝑡) and 𝜔₁(𝑡) are the inputs of the nonlinear static block that represents the cutting process; this subsystem generates the depth of cut 𝑑(𝑡) according to (13). Both the torque-on-bit 𝑇 and weight-on-bit 𝑊 are proportional to 𝑑; they are fed back to rotational motion and vertical motion subsystems, respectively.

Figure 4 

The block diagram of the drilling system.

Our goal is to design a control system that maintains a desired rate of drilling. Specifically, we are looking for the control algorithm for the armature voltage 𝑉 that would guarantee that the velocity of the vertical penetration 𝑣(𝑡) tends asymptotically to an arbitrary positive desired value 𝑣_ref>0. We start designing a control algorithm by considering the equation of vertical motion (15) in some detail.

3.1. Control of the Vertical Motion of a Drill Bit

The vertical motion of the drilling system is described by (15). For convenience, this equation is rewritten below in a slightly modified form, as follows: ̇𝐾𝑣=−^𝑓𝑀𝑊𝑣−^𝑠−𝐻₀𝑀−𝑊𝑀. (16) The idea of the controller developed in this work is to use the weight-on-bit 𝑊 as the control input to the vertical motion subsystem (16). More specifically, combining formulas (7) and (13), one get the following expression for 𝑊: 𝑊=𝑎𝜁𝜖2𝜋𝜔₁𝑣, (17) which essentially indicates that 𝑊 is proportional to the vertical velocity 𝑣(𝑡) and inversely proportional to the angular velocity of the rotational motion 𝜔₁(𝑡). Substituting the last formula into (16), one gets ̇𝑊𝑣=^𝑠−𝐻₀𝑀−1𝑀𝑎𝜁𝜖2𝜋𝜔₁+𝐾_𝑓𝑣. (18) Equation (18) is a linear differential equation with respect to 𝑣 which, assuming 𝜔₁>0, has one stable equilibrium 𝑣=𝑣₀ defined by the formula 𝑊_𝑠−𝐻₀𝑀−1𝑀𝑎𝜁𝜖2𝜋𝜔₁+𝐾_𝑓𝑣₀=0. (19) Solving the above equation with respect to 𝑣₀, one gets 𝑣₀=𝑊_𝑠−𝐻₀𝑎𝜁𝜖2𝜋/𝜔₁+𝐾_𝑓. (20) The above equation (20) indicates that the location of the stable equilibrium 𝑣=𝑣₀ of the vertical motion subsystem (16) can be controlled if one can control the rotational velocity 𝜔₁. Specifically, (20) defines one-to-one correspondence between 𝜔₁ from the range (0,+∞) and 𝑣₀ from the range (0,(𝑊_𝑠−𝐻₀)/𝐾_𝑓). In particular, for any given 𝑣_ref∈(0,(𝑊_𝑠−𝐻₀)/𝐾_𝑓), there exists an unique 𝜔_ref∈(0,+∞) such that if the angular velocity satisfies 𝜔₁(𝑡)≡𝜔_ref, then 𝑣_ref is a globally exponentially stable equilibrium of the translational dynamics (16). For a given 𝑣_ref∈(0,(𝑊_𝑠−𝐻₀)/𝐾_𝑓), the corresponding 𝜔_ref can be found using formula (20), as follows: 𝜔_ref=2𝜋𝑎𝜁𝜖𝑊_𝑠−𝐻₀/𝑣_ref−𝐾_𝑓. (21) Therefore, the control goal of stabilization of the vertical penetration velocity 𝑣(𝑡)→𝑣_ref can be achieved by designing a controller for rotational motion that guarantees a sufficiently fast convergence of 𝜔₁(𝑡)→𝜔_ref. The design of such a controlled is presented in the next section.

3.2. Stabilization of the Angular Velocity of the Drilling System

The rotational dynamics of the drilling system together with the electric drive are described by (5), which is repeated below for convenience, 𝜔̇₁𝜙̇𝜔̇₂𝐼̇=−𝑐₁𝐽₁𝑘𝐽₁000−1010−𝑘𝐽₂−𝑐₂𝐽₂𝐾𝑛𝐽₂00−𝐾𝑛𝐿−𝑅𝐿𝜔₁𝜙𝜔₂𝐼+0001𝐿𝑉+−1𝐽₁000𝑇. (22) The above system has one control input which is the armature voltage of the electric drive 𝑉 and one disturbance input which is the torque-on-bit 𝑇. Our objective in this section is to design a control law for 𝑉 which would track the reference angular velocity of the drill 𝜔₁→𝜔_ref while rejecting the disturbance 𝑇.

To solve the control problem formulated above, one can use the approach to feedforward robust servo control problem presented in [14, 15]. Below, the above approach is described in a simplified manner which, however, serves our purpose well. Consider a linear time invariant system of the form 𝑥̇=𝐴𝑥+𝐵𝑢+𝐷𝑤,𝑦=𝐶𝑥+𝐹𝑢+𝐻𝑤, (23) where 𝑥∈ℝ^𝑛 is the state, 𝑢∈ℝ^𝑚 is the control input, 𝑦∈ℝ^𝑝 is the output, 𝑤∈ℝ^𝑟 are the disturbances, and 𝐴, 𝐵, 𝐶, 𝐷, 𝐹, and 𝐻 are matrices of appropriate dimensions. Consider a control problem described as follows. Suppose the disturbances 𝑤(𝑡) are measurable. Given a desired output signal 𝑦_ref(𝑡), design a control algorithm that guarantees 𝑦(𝑡)→𝑦_ref(𝑡) as 𝑡→+∞. This problem was addressed in [14, 15] in a very general setting. In this work, a simple case is addressed where both 𝑦_ref and 𝑤(𝑡) are assumed to be constant signals, 𝑦_ref(𝑡)≡𝑦_ref and 𝑤(𝑡)≡𝑤_𝑚. In this case, the following two conditions are necessary and sufficient for the existence of a linear time-invariant controller that solves the above described problem.(i)The pair (𝐴,𝐵) is stabilizable, which means that rank𝐵,𝐴𝐵,𝐴²𝐵,…,𝐴^𝑛−1𝐵=𝑛; (24)(ii)Consider rank𝐴𝐵𝐶𝐹=𝑛+𝑝. (25)

If the above two conditions hold (and only in this case), the linear time-invariant controller that solves the above described problem is given according to the formula 𝑢=𝐾𝑥+𝒢^†𝑦_ref+𝒢^∗𝑤_𝑚, (26) where 𝐾∈ℝ^𝑛×𝑛 is the feedback gain matrix which is to be chosen such that 𝐴−𝐵𝐾 is stable and has the required dynamic properties 𝒢=−𝐶(𝐴−𝐵𝐾)⁻¹𝒢𝐵,_∗=𝒢^†𝐶(𝐴−𝐵𝐾)⁻¹𝐷, (27)(28) where 𝒢^† is the Moore-Penrose pseudoinverse of the matrix 𝒢 in (27), defined by the formula 𝒢^†=𝒢^𝑇𝒢𝒢^𝑇−1. (29)

The above described control approach can be applied to the problem of stabilization of the angular velocity of drilling as follows. Equations (22), which describe the rotational dynamics of a drilling system, can be rewritten in the form (23), where 𝑥∶=𝜔₁𝜙𝜔₂𝐼^𝑇∈ℝ⁴, 𝑢∶=𝑉∈ℝ¹, 𝑦∶=𝜔₁∈ℝ¹, 𝑤∶=𝑇∈ℝ¹, and the corresponding matrices are 𝐴∶=−𝑐₁𝐽₁𝑘𝐽₁000−1010−𝑘𝐽₂−𝑐₂𝐽₂𝐾𝑛𝐽₂00−𝐾𝑛𝐿−𝑅𝐿0001,𝐵=𝐿,𝐷=−1𝐽₁000,𝐶=1000,𝐹=0,𝐻=0. (30)(31) Below, we consider the drilling system with specific values of the parameters that are listed in Table 1. With these values, the matrices 𝐴, 𝐵, and 𝐷 become ,000000,𝐴∶=−0.11231.264700−10100−0.2231−0.20050.020400−8640−2𝐵=200,𝐷=−0.0027 (32) while 𝐶, 𝐹, and 𝐻 are given by (31).

Table 1 

Numerical values for drilling system parameters.

For the above system, the necessary and sufficient conditions for stabilization (24), (25) are satisfied. Indeed, the stabilizability condition (24) is satisfied since rank𝐵,𝐴𝐵,𝐴²𝐵,…𝐴^𝑛−1𝐵=rank0005.15427300.0000004.075472−8.96804.075−8.968−700.339200−400−34412146307=4. (33) On the other hand, the rank condition (25) is also satisfied because rank𝐴𝐵𝐶𝐹=rank−0.1122991.264706000−101000−0.223113−0.2004720.020377000−8640−220010000=5. (34) Therefore, a controller of the form (26), (27), (28), and (29) guarantees that the angular velocity of the drill approach the reference angular velocity 𝜔₁→𝜔_ref as 𝑡→∞, while rejecting the disturbance 𝑇_𝑏.

The design of controller (26), (27), (28), and (29) begins by choosing the desired location of the closed-loop system's poles. For the purpose of simulations presented below, we consider two specific set of poles. The first set, denoted by 𝑃₁, is chosen as follows: 𝑃₁∶=[−10−2+2𝑖−2−2𝑖−4]. (35) The set 𝑃₁ consists of two real poles and two complex conjugate poles. On the other hand, the set 𝑃₂ contains only poles on the real axis as follows: 𝑃₂=[−5.5−2−4.5−1]. (36) The feedback gain matrix 𝐾₁ such that the poles of 𝐴−𝐵𝐾₁ are located according to 𝑃₁ is 𝐾₁=[32.2457.45−19.410.0784]. (37) The coefficients 𝒢^∗, 𝒢^† in (26) are calculated according to the formulas (27)–(29); the results are 𝒢^∗₁=0.123497,𝒢^†₁=60.0844. (38) On the other hand, the feedback matrix 𝐾₂ such that the poles of 𝐴−𝐵𝐾₂ are located according to 𝑃₂ is 𝐾₂=[−5.16716.943−30.620.0534]. (39) The corresponding coefficients 𝒢^∗₂, 𝒢^†₂ are 𝒢^∗₂=0.037286,𝒢^†₂=9.603682. (40)

3.3. Rock Stiffness Estimation

In the controller design presented above, it was assumed that the “hardness” of the rock, represented by the intrinsic specific energy 𝜖, is constant and exactly known. This knowledge of 𝜖 was used explicitly in the controller design, in particular, in formula (21). In practical geological drilling, however, the hardness of different layers of rock lying underneath the surface can be different and usually is not exactly known beforehand. More specifically, different characteristics of the rock, such as hardness, density and porosity, typically remain constant through each layer, but differs from layer to layer. On the other hand, control engineers frequently deal with the problem of designing a controller without a priori knowledge of the exact values of one or more parameters involved in the process. Often, the processes can be robustly controlled without the actual knowledge of some of the parameters. In other cases, the unknown parameters can be identified using specially designed estimators. Below, a simple online estimator of the rock intrinsic specific energy 𝜖 is designed following the methods described in [16], and the resulting estimate is then used in the controller for for drilling system.

Specifically, during the cutting process, the torque-on-bit 𝑇 is produced by bit rock interaction, according to the formula 1𝑇=2𝑎²𝜖𝑑, (41) where 𝑎 is the radius of drill bit, 𝑑 is the depth of cut, and 𝜖>0 is the intrinsic specific energy. The intrinsic specific energy 𝜖>0 depends on the properties of the media and typically unknown beforehand. However, since the torque on bit 𝑇(𝑡) can typically be measured with advanced transducers located in the bottom hole assembly [17], 𝑎>0 is constant and known, and 𝑑(𝑡) can be calculated according to formula (13), one can use the method described in the previous section to design an online estimation scheme for 𝜖. In particular, considering (1/2)𝑎²𝑑(𝑡) as the input and torque-on-bit 𝑇 as the measured output, one can follow the procedure described in the previous section to design an estimator for an unknown parameter 𝜖. The predicted torque-on-bit 𝑇_𝑏 is defined according to the formula 𝑇(𝑡)1∶=2𝑎²𝜖𝑑̂(𝑡), (42) where 𝜖(̂𝑡) is the current estimate of actual rock strength 𝜖. The algorithm for online estimation of the intrinsic specific energy 𝜖 has a form 𝜖̇̂=𝛾₀𝑇−𝑇12𝑎²𝑑, (43) where 𝛾₀>0 is an arbitrary gain.

A natural question regarding the algorithm (43) is if it guarantees the convergence of the parameter estimate to the true value of the parameter 𝜖; mathematically, is 𝜖(̂𝑡)→𝜖 as 𝑡→+∞. It is known [16], that the convergence can be guaranteed if the “input” signal (1/2)𝑎²𝑑(𝑡) is persistently exciting. A signal 𝑢(𝑡) is said to be persistently exciting which is to say that there exist 𝛼₀>0, 𝑇₀>0 such that the inequality ^𝑡+𝑇0_𝑡𝑢²(𝜏)𝑑𝜏≥𝛼₀𝑇₀ (44) holds for all 𝑡. In particular, 𝑢(𝑡) is persistently exciting if 𝑢²(𝑡)≥𝛼₀ for all 𝑡. Since 𝑑(𝑡) is the depth of cut, we see that, during normal cutting process, 𝑑(𝑡)≥𝑑₀>0, which results in persistent excitation of the input (1/2)𝑎²𝑑(𝑡). The parameter convergence 𝜖(̂𝑡)→𝜖, therefore, is guaranteed during normal cutting process. This is also confirmed by the simulation results presented below.

The obtained estimate of the rock strength 𝜖̂ is then used in the control algorithm. Specifically, in the original formulation of the control algorithm, for a given reference vertical velocity 𝑣_ref, the reference rotational velocity 𝜔_ref is calculated according to formula (21), which depends on the parameter 𝜖. In case 𝜖 is unknown, it is substituted by its estimate 𝜖(̂𝑡) obtained above. The new formula for 𝜔_ref has the form 𝜔_ref𝜖̂∶=2𝜋𝑎𝜁𝑊_𝑠−𝐻₀/𝑣_ref−𝐾_𝑓. (45) The obtained estimate of the rock stiffness 𝜖̂ will also be used to update the stiffness of the virtual spring in the haptic teleoperator drilling system described below.

3.4. Simulation Results

In this Section, some results of simulations of the drilling control system with intrinsic specific energy estimator are presented. The vertical motion of the drilling system is described by (15), and it is interconnected with the rotational dynamics (5) through nonlinear equation (13) that describes the depth of cut 𝑑(𝑡). For a given reference velocity of the vertical penetration 𝑣_ref>0, the corresponding reference rotational velocity 𝜔_ref is calculated according to formula (45). The controller (21), (26)–(29) has been implemented to guarantee that the angular velocity of the drill bits 𝜔₁(𝑡) tracks 𝜔_ref, which in turn stabilizes the vertical penetration velocity 𝑣(𝑡) converges to 𝑣_ref. The algorithm (43) provides an estimate of the intrinsic specific energy parameter 𝜖̂ which is then used in the calculation of the reference angular velocity according to formula (45). Specific values of the parameters appearing in these equations are given in Table 1. The simulations are carried out using MATLAB, where the integration step for each simulation is equal to 0.005 s. The feedback gain matrix is chosen 𝐾=𝐾₁, where 𝐾₁ is defined by (37).

In the simulations described below, the performance of the system was evaluated for different values of actual intrinsic specific energy 𝜖, different gains 𝛾₀ and different values of the applied weight 𝑊₀∶=𝑊_𝑠−𝐻₀. Figures 5 and 6 show the response of the vertical penetration velocity 𝑣(𝑡), the intrinsic specific energy estimate 𝜖(̂𝑡), the torque-on-bit 𝑇(𝑡), the predicted value of the torque-on-bit 𝑇(𝑡), and the rotational velocity 𝜔₁(𝑡), all for the case where the applied weight on bit 𝑊₀=5000 N, the intrinsic specific energy 𝜖=20 MPa, and the desired vertical velocity 𝑣_ref is set to 0.005 m/s. The estimator gain is set to 𝛾₀=5⋅10⁹. The plots show that 𝑣(𝑡) converges to 𝑣_ref in less than 8 sec whereas the estimate 𝜖(̂𝑡) converges to the actual value of 𝜖 in less than 4 sec. Figures 7 and 8 show the output responses of described parameters where 𝑊₀=2500 N and the desired vertical velocity 𝑣_ref is set to 0.01 m/s. It can be clearly seen that the convergence becomes slower with reducing the applied weight on the drill string 𝑊₀; specifically, both 𝑣(𝑡) and 𝜖(̂𝑡) approach their reference values in about 12 sec. Reducing 𝑊₀ also results in that 𝜔_1ref increases, the steady-state value of 𝑇_𝑏(𝑡) drops to around 200 N, and the steady-state value of 𝑑(𝑡) also drops to less than 2 mm. On the other hand, Figures 9 and 10 demonstrate the response of the system with the same parameters except the intrinsic specific energy 𝜖 is reduced to 5 MPa. This results in decreased convergence time for 𝑣(𝑡) and 𝜖(̂𝑡). The steady state value of rotational velocity 𝜔₁(𝑡) is also decreased to under 10 rad/s, and steady state value of the depth of cut 𝑑(𝑡) is increased to 6.5 mm. Figures 11 and 12 present the output response for the case where the estimator gain is decreased to 𝛾₀=5⋅10⁸, while the rest of the parameters are the same as in the last simulation except the intrinsic specific energy is set to 𝜖=10 MPa. The resulting response is predictably characterized by much slower convergence, which takes about 25 sec for 𝑣(𝑡) and 𝜖(̂𝑡) to approach their steady-state values. Finally, Figures 13 and 14 correspond to to the case where 𝑊₀=5000 N, 𝜖=20 MPa, and 𝛾₀=1⋅10⁸.

(a)

(b)

Figure 5 

Response of the vertical velocity 𝑣(𝑡) (a) and the intrinsic specific energy estimate 𝜖(̂𝑡) (b) for 𝑊0=5000 N, 𝜖=20 MPa, and 𝛾0=5⋅109.

(a)

(b)

(c)

Figure 6 

Response of rotational velocity 𝜔1(𝑡) (a), torque-on-bit 𝑇(𝑡) vresus estimated torque-on-bit 𝑇(𝑡) (b), and the depth of cut 𝑑(𝑡) (c) for 𝑊0=5000 N, 𝜖=20 MPa, and 𝛾0=5⋅109.

(a)

(b)

Figure 7 

Response of the vertical velocity 𝑣(𝑡) (a) and the intrinsic specific energy estimate 𝜖(̂𝑡) (b) for 𝑊0=2500 N, 𝜖=20 MPa, and 𝛾0=5⋅109.

(a)

(b)

(c)

Figure 8 

Response of rotational velocity 𝜔1(𝑡) (a), torque-on-bit 𝑇(𝑡) versus estimated torque-on-bit 𝑇(𝑡) (b), and the depth of cut 𝑑(𝑡) (c) for 𝑊0=2500 N, 𝜖=20 MPa, and 𝛾0=5⋅109.

(a)

(b)

Figure 9 

Response of the vertical velocity 𝑣(𝑡) (a) and the intrinsic specific energy estimate 𝜖(̂𝑡) (b) for 𝑊0=2500 N, 𝜖=5 MPa, and 𝛾0=5⋅109.

(a)

(b)

(c)

Figure 10 

Response of rotational velocity 𝜔1(𝑡) (a), torque-on-bit 𝑇(𝑡) versus estimated torque-on-bit 𝑇(𝑡) (b), and the depth of cut 𝑑(𝑡) (c) for 𝑊0=2500 N, 𝜖=5 MPa, and 𝛾0=5⋅109.

(a)

(b)

Figure 11 

Response of the vertical velocity 𝑣(𝑡) (a) and the intrinsic specific energy estimate 𝜖(̂𝑡) (b) for 𝑊0=2500 N, 𝜖=10 MPa, and 𝛾0=5⋅108.

(a)

(b)

(c)

Figure 12 

Response of rotational velocity 𝜔1(𝑡) (a), torque-on-bit 𝑇(𝑡) versus estimated torque-on-bit 𝑇(𝑡) (b), and the depth of cut 𝑑(𝑡) (c) for 𝑊0=2500 N, 𝜖=10 MPa, and 𝛾0=5⋅108.

(a)

(b)

Figure 13 

Response of the vertical velocity 𝑣(𝑡) (a) and the intrinsic specific energy estimate 𝜖(̂𝑡) (b) for 𝑊0=2500 N, 𝜖=20 MPa, and 𝛾0=108.

(a)

(b)

(c)

Figure 14 

Response of rotational velocity 𝜔1(𝑡) (a), torque-on-bit 𝑇(𝑡) versus estimated torque-on-bit 𝑇(𝑡) (b), and the depth of cut 𝑑(𝑡) (c) for 𝑊0=2500 N, 𝜖=20 MPa, and 𝛾0=108.

Overall, simulation results show that the control system with intrinsic specific energy estimation demonstrate good stability and performance characteristics for a wide range of the parameters. In particular, the vertical velocity converges to the desired value, and the estimate of the intrinsic specific energy 𝜖(̂𝑡) converges to an actual value of 𝜖.

4. Telerobotic Drilling System with Haptic Feedback

In this section, a telerobotic drilling system with haptic feedback is designed and experimentally evaluated. Haptics can be defined as the physical or virtual interaction through touch sensation for the purpose of perception and manipulation of objects [18, 19]. Haptic feedback provides the operator with kinaesthetic clues of the physical features of virtual or real remote environment. The structure of a telerobotic drilling system with haptic feedback is shown in Figure 15. In this system, the human operator controls the drilling process using a haptic device. Specifically, the position of an end-effector of the haptic device defines the reference vertical velocity of the drilling. The reference vertical velocity is then transmitted to the drilling control system, designed in above in Section 3, which stabilizes the actual vertical penetration velocity to the level equal to the reference vertical velocity. On the other hand, an estimate of the intrinsic specific energy 𝜖(𝑡), which is generated online by an estimator described in above in Section 3.3, is sent back to the haptic device. The end-effector of the haptic device interacts with a virtual spring of variable stiffness; the stiffness of this virtual spring is updated in real time proportionally to the current estimate of the intrinsic specific energy 𝜖(𝑡). Thus, the telerobotic drilling system provides haptic feedback to the human operator which creates an intuitive feeling of the hardness of the remotely drilled material.

Figure 15 

The structure of a telerobotic drilling system.

4.1. Experimental Setup

The above described telerobotic drilling system is implemented in a semiexperimental setup, as follows. The setup consists of a PC based on Intel Pentium 4 processor with operational frequency of 1 GHz and RAM of 1 GB, and a PHANTOM Omni Haptic device, a product from SensAble Technologies Inc. The PHANTOM Omni Haptic device is designed for kinematic interaction with the virtual or real environment while providing the kinesthetic feedback to the operator. The device is equipped with a pen-based stylus, and is able to provide three degrees-of-freedom force feedback. The human operator uses the haptic device to (i) generate a desired vertical velocity 𝑣_ref(𝑡) which is used as an input to the drilling control system, and (ii) to haptically perceive the stiffness 𝜖 of the rock layers. The remaining parts of the above described telerobotic system, including the drill string and drive system, the drilling process, as well as the control and estimation algorithms, are simulated in real time in virtual environment which is implemented using the Open Haptics tool kit and Microsoft Visual C++.

The human operator assigns the desired velocity 𝑣_ref(𝑡) by controlling the position of the end-effector of the PHANTOM device along its vertical (𝑌) axis. More exactly, a specific range along the 𝑦-axis is assigned to each desired vertical velocity, as follows: 𝑣_ref(𝑡)=0.001 m/s if the position of stylus 𝑦_𝑛(𝑡) is ≥80 mm; 𝑣_ref(𝑡)=0.003 m/s if the position of stylus 𝑦_𝑛(𝑡) is ≤80 mm and ≥60 mm;  𝑣_ref(𝑡)=0.005 m/s if the position of stylus 𝑦_𝑛(𝑡) is ≤60 mm and ≥40 mm;  𝑣_ref(𝑡)=0.008 m/s if the position of stylus 𝑦_𝑛(𝑡) is ≤40 mm and ≥25 mm;  𝑣_ref(𝑡)=0.01 m/s if the position of stylus 𝑦_𝑛(𝑡) is ≤25 mm and ≥10 mm;  𝑣_ref(𝑡)=0.015 m/s if the position of stylus 𝑦_𝑛(𝑡) is ≤10 mm and ≥0 mm;  𝑣_ref(𝑡)=0.018 m/s if the position of stylus 𝑦_𝑛(𝑡) is ≤0 mm.

Another function of the haptic device is to allow the human operator to feel the stiffness of the rocks. As explained above, this is achieved by updating the stiffness of the virtual spring proportionally to the current estimate of the rock stiffness (intrinsic specific energy 𝜖(̂𝑡)). The coefficient of proportionality between the estimate of the intrinsic specific energy (with units of Pascals) and the stiffness of the virtual spring (with units of is N/m) is set in our experiments equal to 10⁻⁷. The feedback force 𝐹_est(𝑡) due to the virtual spring is therefore calculated according to the formula 𝐹_est(𝑡)=10⁻⁷⋅𝜖̂(𝑡)⋅𝑦_𝑛(𝑡). (46)

4.2. Experimental Results

In this Section, some experimental results are presented. In these experiments, we have attempted to simulate a real drilling case scenario where the composition, characteristics, and types (which all contribute to the intrinsic specific energy) of various rock strata vary at different depths during drilling. Specifically, in every experiment, several layers of rocks with different intrinsic specific energy 𝜖 ranging from 4 MPa to 60 MPa are simulated. In all experiments presented here, the applied weight 𝑊₀=5000 N; the rest of the parameters, if not explicitly mentioned, are same as in Section 3.4.

In the experiment shown in Figures 16, 17, and 18, three layers of rocks with different intrinsic specific energy 𝜖 are simulated. The top layer has the stiffness of 5 MPa and its thickness is 20 cm from the surface. The second layer has a stiffness value of 12 MPa and lies between 20 cm and 30 cm from the surface (total thickness is 10 cm). The third layer starts at the depth of 30 cm and continues downward. It has a stiffness value of 20 MPa. The experiment is performed with the estimator gain 𝛾₀=10⁹. Figure 16 shows the actual intrinsic specific energy 𝜖(𝑡) and its estimate 𝜖(̂𝑡) on the top graph, and the reflected force 𝐹_est(𝑡) on the bottom graph. Due to high estimator gain, 𝜖(̂𝑡) quickly tracks 𝜖(𝑡) for all three layers as the drill bit progressed cutting through these layers. Figure 17 shows the vertical velocity 𝑣_out(𝑡) and the reference vertical velocity 𝑣_ref(𝑡) on the top graph, and the reference rotational velocity 𝜔_1𝑑(𝑡) and the actual drill bit rotational velocity 𝜔₁(𝑡) at the bottom graph. Finally, Figure 18 shows the behaviour of the actual torque-on-bit 𝑇(𝑡) and the estimated torque 𝑇(𝑡), along with depth of cut 𝑑(𝑡). These plots show that the system is stable and demonstrates good performance; in particular, all the output variables track their desired (reference) trajectories.

(a)

(b)

Figure 16 

Experiment 1: Actual stiffness 𝜖(𝑡) versus the estimated stiffness 𝜖(̂𝑡) (a); the reflected force 𝐹est(𝑡) (b).

(a)

(b)

Figure 17 

Experiment 1: Output vertical velocity 𝑣out(𝑡) versus reference vertical velocity 𝑣ref(𝑡) (a); output rotational velocity of the drill bit 𝜔1(𝑡) versus reference rotational velocity 𝜔1𝑑(𝑡) (b).

(a)

(b)

Figure 18 

Experiment 1: Torque-on-bit 𝑇(𝑡) versus estimated torque-on-bit 𝑇(𝑡) (a); depth of cut 𝑑cut(𝑡) (b).

Another set of experimental results is presented in Figures 19–21, where the estimator gain is increased to 𝛾₀=5⋅10⁹, and the depth of the rock layers and their corresponding stiffness values have been altered. Specifically, the first rock layer has depth 20 cm and the intrinsic specific energy 𝜖 is set to 20 MPa for this layer. Second layer lies between 20 cm and 40cm with 𝜖=40 MPa. The third layer lies below 40 cm, and its intrinsic specific energy 𝜖=60 MPa. Figure 19 shows the corresponding plots of 𝜖(𝑡), 𝜖(̂𝑡) and 𝐹_est(𝑡). Figure 20 shows the response of 𝑣_ref(𝑡) and 𝑣_out(𝑡) on the top graph, and the responses of 𝜔_1𝑑(𝑡) and 𝜔₁(𝑡) on the bottom graph, respectively. The response of 𝑇(𝑡) and 𝑇(𝑡) along with 𝑑(𝑡) are shown in Figure 21. Overall, our experiments demonstrate stability and good performance of the designed telerobotic drilling system with haptic feedback, for a wide range of parameters and control gains.

(a)

(b)

Figure 19 

Experiment 2: Actual stiffness 𝜖(𝑡) versus estimated stiffness 𝜖(̂𝑡) (a); reflected force 𝐹est(𝑡) (b).

(a)

(b)

Figure 20 

Experiment 2: Output vertical velocity 𝑣out(𝑡) versus reference vertical velocity 𝑣ref(𝑡) (a); output rotational velocity of the drill bit 𝜔1(𝑡) versus reference rotational velocity 𝜔1𝑑(𝑡) (b).

(a)

(b)

Figure 21 

Experiment 2: Torque-on-bit 𝑇(𝑡) versus estimated torque-on-bit 𝑇(𝑡) (a); depth of cut 𝑑cut(𝑡) (b).

5. Conclusions

This paper deals with control design for a teleoperator system with haptic feedback for an oil well drilling process. A mathematical model of the drilling process was described, and the control algorithm was designed that guarantees the convergence of the vertical penetration velocity to an arbitrary reference value. The control algorithm has a cascaded structure, where the velocity of vertical penetration is controlled indirectly through stabilization of the rotational motion of the drill bit. In order to guarantee the convergence of the angular velocity to a desired value in the presence of disturbances in the form of torque-on-bit, a robust servo controller was designed. However, the design of such controller depends on the parameter of environment called the intrinsic specific energy, which is generally unknown beforehand. To solve this issue, an online parameter estimator was designed that provides an estimate of the intrinsic specific energy. This estimate is substituted for the actual value of the parameter in the control algorithm, and the corresponding adaptive control system is evaluated through simulations. Finally, a telerobotic drilling system with haptic feedback is designed and verified through semi-experiments. The haptic feedback for the human operator is provided by creating a virtual spring that interacts with the haptic device; the stiffness of the spring is adjusted in real time depending on the current estimate of the intrinsic specific energy. Semi-experiments are conducted using PHANTOM Omni Haptic device, where the drilling process model is implemented in C++ environment, and the haptic feedback is provided to the human operator.

There exists a number of challenges associated with the real-life drilling operation that were not addressed in our paper. In particular, the frictional forces at the contact were neglected in our analysis, while in reality they may play significant role in the drilling process. Also, in real-life drilling systems, the rotational velocity and the penetration rate are typically measured at the surface while the torque-on-bit should ideally be measured above the bit; thus, there exists a problem of synchronizing the data obtained at the surface with those obtained at the bit. The issue of communication delay between the haptic device and the drilling process is also not addressed. These, as well as more detailed experimental evaluation of the designed system, are the topics for future research.

Conflict of Interests

The authors of the paper do not have a direct financial relation with the commercial identity mentioned in the paper that might lead to a conflict of interests for any of the authors.

Acknowledgment

This work was supported by the Natural Sciences and Engineering Research Council (NSERC) of Canada under Discovery Grant RGPIN1510.