Method of and a device for estimating down hole speed and down hole torque of borehole drilling equipment while drilling, borehole equipment and a computer program product

Abstract

A method of and a device, borehole drilling equipment and a computer program product for estimating at least one of down hole speed and down hole torque of borehole drilling equipment while drilling a borehole in an earth formation. The borehole equipment comprising a rotational drive system, a drill string, a bottom hole assembly comprising a drill bit, a top end coupled to the rotational drive system, and a speed controller arranged for controlling rotational top end drive speed by commanding a top end drive torque. The method being performed in a computer based on a two-port spectral domain transfer matrix computational model of the borehole drilling equipment. While drilling, rotational top end drive speed and top end drive torque are obtained and down hole speed and down hole torque of the borehole drilling equipment are calculated by reducing spectral bandwidth of and applying a time delay term to the calculations enabling a causal time domain solution.

Claims

1. A method of estimating at least one of down hole speed and down hole torque of borehole drilling equipment while drilling a borehole in an earth formation, said borehole drilling equipment comprising a rotational drive system, a drill string having a bottom hole assembly comprising a drill bit and a top end coupled to said rotational drive system, and a speed controller arranged for controlling rotational top end drive speed by commanding a top end drive torque, the method being performed in a computer and comprising the steps of: obtaining transfer matrix components of a two-port spectral domain transfer matrix computational model of said borehole drilling equipment, said computational model including damping properties of said borehole drilling equipment, wherein said down hole speed and down hole torque being dependent variables and said rotational top end drive speed and top end drive torque being independent variables; obtaining, while drilling, said rotational top end drive speed and top end drive torque, and calculating, from said obtained transfer matrix components, said top end drive speed and said top end drive torque at least one of said down hole speed and down hole torque of said borehole drilling equipment, wherein said step of calculating comprises: reducing spectral bandwidth of said calculation, said reduced bandwidth being selected to include a limited number of higher order spectral modes in said calculation, applying, in said spectral domain, a time delay term to said calculation, said time delay term being selected to enable a causal time domain solution, and calculating at least one of said down hole speed and down hole torque of said borehole drilling equipment from applying said reduced bandwidth and added time delay term.

2. The method according to claim 1, wherein said spectral bandwidth reduction and said time delay term are applied by modifying said transfer matrix components by: reducing spectral bandwidth of said transfer matrix components, said reduced bandwidth being selected to include a limited number of higher order spectral modes in said transfer matrix components, adding, in said spectral domain, a time delay term to each of said transfer matrix components, said time delay being selected to enable a causal time domain solution, and calculating at least one of said down hole speed and down hole torque of said borehole drilling equipment using said modified transfer matrix components.

3. The method according to claim 2, wherein at least one of said down hole speed and down hole torque of said borehole drilling equipment are calculated, by said computer, in a causal convolution operation processing said obtained rotational top end drive speed and top end drive torque and impulse response functions of said modified transfer matrix components over a period of time including one round trip time in said drill string.

4. The method according to claim 2, wherein at least one of said down hole speed and down hole torque are provided, by said computer, for a plurality of different positions along said drill string length or borehole by calculating, for each position, impulse response functions of said modified transfer matrix transfer function components, and by processing in a causal convolution operation from each of said calculated impulse response functions for a respective position and said obtained rotational top end drive speed and top end drive torque at least one of said down hole speed and down hole torque at a respective position.

5. The method according to claim 4, wherein said at least one of said down hole speed and down hole torque provided for said plurality of positions are represented graphically, while drilling, at a display device in a spatial diagram as a function of time and position.

6. The method according to claim 1, wherein said spectral bandwidth is reduced based on at least one of damping properties and length of said drill string.

7. The method according to claim 1, wherein said spectral bandwidth is reduced by applying an nth order low-pass Bessel filter operation, wherein n≥4, and tapering, wherein said time delay term further comprises a time delay part representative of a time delay caused by said Bessel filter operation and n is a positive integer.

8. The method according to claim 1, wherein said time delay term equals at least 1 second per 3000 meter of drill string length.

9. The method according to claim 1, wherein said down hole speed and down hole torque are rotational bit speed and bit torque, respectively.

10. The method according to claim 1, wherein said two-port spectral domain transfer matrix computational model of said borehole drilling equipment comprises a two-port spectral domain transfer matrix, composed of a single or composite symmetrical spectral domain transmission line equivalent drill string computational model including drill string inertia, stiffness and damping properties of said drill string in said borehole, bottom hole assembly inertia and top end drive inertia.

11. The method according to claim 10, wherein said two-port spectral domain transfer matrix computational model is represented by: $\left(\begin{array}{c}{\omega }_{\mathrm{dh}}\\ T_{\mathrm{dh}}\end{array}\right)=\left[\begin{array}{cc}{M}_{11}\left(s\right)& M_{12}\left(s\right)\\ M_{21}\left(s\right)& M_{22}\left(s\right)\end{array}\right]\left(\begin{array}{c}{\omega }_{\mathrm{td}}\\ T_{\mathrm{ref}}.Math.e^{-s{\tau }_d}\end{array}\right)$ wherein: M11(s), M12(s), M21(s) and M22(s) represent bandwidth limited and time delay term modified spectral transfer matrix transfer function components, ω.sub.dh represents down hole speed, T.sub.dh represents down hole torque, ω.sub.td represents obtained top end drive speed, T.sub.ref represents commanded top end drive torque, τ.sub.d represents speed controller time delay, and s represents Laplace transform operator.

12. The method according claim 1, wherein the method is performed in a computer separate from said borehole drilling equipment.

13. A computer readable product, storing computer program code on a non-transitory medium that can be read by a computer, wherein said computer program code is arranged to carry out the method according to claim 1 when said computer is operatively connected to borehole drilling equipment and said program code is executed by said computer.

14. A device for estimating at least one of down hole speed and down hole torque of borehole drilling equipment for drilling a borehole in an earth formation, said borehole drilling equipment comprising a rotational drive system, a drill string having a bottom hole assembly comprising a drill bit and a top end coupled to said rotational drive system, and a speed controller arranged for controlling rotational top end drive speed by commanding a top end drive torque, said device comprising a computer controlled drilling operation observer system arranged for: obtaining transfer matrix components of a two-port spectral domain transfer matrix computational model of said borehole drilling equipment, said computational model including damping properties of said borehole drilling equipment, and wherein said down hole speed and down hole torque being dependent variables and said rotational top end drive speed and top end drive torque being independent variables; obtaining, while drilling, said rotational top end drive speed and top end drive torque, and calculating, from said obtained transfer matrix components, said top end drive speed and said top end drive torque at least one of said down hole speed and down hole torque of said borehole drilling equipment, wherein said calculating comprises: reducing spectral bandwidth of said calculation, said reduced bandwidth being selected to include a limited number of higher order spectral modes in said calculation, applying, in said spectral domain, a time delay term to said calculation, said time delay term being selected to enable a causal time domain solution, and calculating at least one of said down hole speed and down hole torque of said borehole drilling equipment from applying said reduced bandwidth and added time delay term.

15. The device according to claim 14, wherein said computer controlled drilling operation observer system is arranged for applying said spectral bandwidth reduction and said time delay term by modifying said transfer matrix components by: reducing spectral bandwidth of said transfer matrix components, said reduced bandwidth being selected to include a limited number of higher order spectral modes in said transfer matrix components; adding, in said spectral domain, a time delay term to each of said transfer matrix components, said time delay being selected to enable a causal time domain solution, and calculating at least one of said down hole speed and down hole torque of said borehole drilling equipment using said modified transfer matrix components.

16. The device according to claim 15, wherein said computer controlled drilling operation observer system is arranged for calculating, for each position of a plurality of different positions along said drill string length or borehole, impulse response functions of said modified transfer matrix transfer function components, and by processing in a causal convolution operation from each of said calculated impulse response functions for a respective position and said obtained rotational top end drive speed and top end drive torque at least one of said down hole speed and down hole torque at a respective position.

17. The device according to claim 16, further comprising a display device, wherein said computer controlled drilling operation observer system is arranged for graphically representing said at least one of said down hole speed and down hole torque provided for said plurality of positions, while drilling, at said display device in a spatial diagram as a function of time and position, in particular in a speed wave pattern or torque wave pattern diagram.

18. The device according to claim 14, wherein said computer controlled drilling operation observer system is arranged for reducing said spectral bandwidth by applying an nth order low-pass Bessel filter operation, wherein n≥4, and tapering, and adding to said time delay term a time delay part representative of a time delay caused by said Bessel filter operation, and wherein n is a positive integer.

19. The device according to claim 14, operatively connected to or incorporated in said speed controller for obtaining said top end drive speed and top end drive torque.

20. Borehole drilling equipment for drilling a borehole in an earth formation, said borehole drilling equipment comprising a rotational drive system, a drill string having a bottom hole assembly comprising a drill bit and a top end coupled to said rotational drive system, a speed controller arranged for controlling a commanded torque and rotational drive speed of said top end, and a device for estimating at least one of down hole speed and down hole torque of said borehole drilling equipment in accordance with claim 14.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) FIG. 1 is a very schematic representation of prior art borehole drilling equipment for drilling a borehole in an earth formation.

(2) FIG. 2 shows a system model in the spectral domain of the borehole drilling equipment shown in FIG. 1, with a lumped drill string model.

(3) FIG. 3 shows a distributed drill string model in the spectral domain.

(4) FIG. 4 shows an electrical equivalent circuit diagram of a string fraction of the distributed string model of FIG. 3.

(5) FIG. 5 shows an electrical equivalent transmission line circuit diagram applicable to model a drill string.

(6) FIGS. 6a, 6b, 6c and 6d show graphs of impulse response functions calculated for a practical borehole drilling equipment simulated at a simulation test system of applicant, with the drill string represented by a single segment or section model and a four segment or section model, in accordance with the present disclosure.

(7) FIGS. 7a, 7b and 7c illustrate graphically, over a selected time period, the bit speed and the contributions to the bit speed calculated for the borehole drilling equipment corresponding to FIGS. 6a-6d, based on the four segment drill string impulse response functions and the obtained actual values of the top end drive speed and the commanded top end drive torque at the simulation test system.

(8) FIG. 8 illustrates graphically, on an enlarged time scale with respect to FIGS. 7a-7c, actual measured bit speed values and the bit speed calculated in accordance with the present disclosure from a four segment drill string model and a single segment model of the drill string for the borehole drilling equipment simulated at the simulation test system corresponding to FIGS. 6a-6d.

(9) FIGS. 9a, 9b, 9c and 9d illustrate graphically, over a selected time period, the contributions to the bit torque calculated for the borehole drilling equipment based on the four drill string impulse response functions and the obtained actual values of the top end drive speed and the commanded top end drive torque, and the calculated and measured bit torque simulated at the test system corresponding to FIGS. 6a-6d.

(10) FIGS. 10a and 10b illustrate graphically, at different time instances, the down hole speed without and with applying stick-slip mitigating control, respectively, estimated in accordance with the present disclosure.

(11) FIG. 11 illustrates graphically, in a spatial or 3D diagram, the down hole speed while drilling, estimated in accordance with the present disclosure, for a selected time period before and after applying stick-slip mitigating control.

(12) FIG. 12 illustrates graphically, in a spatial or 3D diagram, the down hole torque, estimated in accordance with the present disclosure, for a selected time period before and after applying stick-slip mitigating control.

(13) FIG. 13 shows a simplified flow chart diagram for automatically determining the down hole speed and/or down hole torque from a computational model of borehole drilling equipment in accordance with the present disclosure.

(14) FIG. 14 is a schematic representation of borehole drilling equipment equipped for operating in accordance with the present disclosure.

DETAILED DESCRIPTION

(15) FIG. 1 shows, in a very schematic manner, a typical borehole drilling equipment 10 of a drilling rig for drilling a borehole in an earth formation. A cutting tool or drill bit 17 connects to a bottom hole assembly, BHA, 11 at a bottom end 13 of a drill string 12. At a top end 14 thereof, the drill string 12 is coupled to a rotational drive system 15, also called top drive or rotary table, which, in turn, is fixed to the surface of an earth formation in which a borehole is to be drilled by a derrick or a chassis of a drilling rig or may be mounted at a ship or the like (not shown).

(16) The drill string 12 comprises lengths of hollow tubulars or drill pipes, threaded together end by end. A typical drill string is several kilometres long, such as 0-10 km, and the drill pipe may have an outer diameter of about 100-300 mm and a wall thickness of about 10-50 mm. The BHA 11 consists of heavier pipes that may have an outer diameter of about 250-500 mm and a wall thickness of about 100 mm, for example. The length of the BHA is typically in the range of 100-300 m. The drill string 12 is very slender compared to its length.

(17) Although not shown, in an actual drilling operation, drilling fluid is pumped through the drill pipes of the drill string 12 towards the drill bit 17 for cooling and lubrication of the drill bit 17. Cuttings from the drilling operation are returned back up to the surface by the drilling fluid flowing through an annulus formed between the outer circumference of the drill string 12 and the borehole (not shown).

(18) The bottom hole assembly 11 comprises several sensors and transmitters 16 and a directional tool (not shown) for directing the bottom hole assembly 11 to drill a borehole in a certain direction in the earth formation, such as vertical, horizontal or deviated at an angle and, of course, combinations thereof.

(19) Drilling data and information are displayed at a console 19 comprising a display or other data output device (not shown) and an input device such as a keyboard, touch screen and the like (not shown) by which, through an intermediate speed controller 20, a driller may control rotational speed of the drive system 15 by inputting tuning parameters for the speed controller 20 and/or setting a torque limit for the drive system 15, for controlling the rotational speed of the drill bit 17.

(20) The drive system 15 comprises a rotary drive system motor 18 to rotate the drill string 12, the BHA 11 and thereby the drill bit 17. Between the drive system motor 18 and the drill string 12 optionally a gearbox 22 may connect, having a particular gear reduction or a range of gear reductions.

(21) Nowadays the drive system motor 18 generally is an electric motor, for example an 800 kW induction motor powered by a power converter. However, the present disclosure is equally applicable to a synchronous machine, a brushed DC machine, diesel engine, a hydraulic motor, or the like. The rotational speed of the drill string 12 may be measured at its top end 14 by a speed indicator or speed sensor 21, for example a speed sensor at the fast rotating motor shaft, the measurement signal of which is input to the speed controller 20. However, in the case of an electric drive system motor, the top end drive speed of the drill string may be obtained, for or with the speed controller 20, in a sensorless manner from the voltage and current supplied to the electric drive system motor 18, for example.

(22) In use, at its top end 14, the drill string 12 can be pulled upwards by hoisting machinery, also called draw-works (not shown). On the bottom end 13, when on-bottom, the BHA 11 is resting at or touches the earth formation by the drill bit 17. This contrary to an off-bottom BHA 11 position, in which the drill bit 17 does not touch the earth formation at the bottom of the borehole. The slender drill pipes of the drill string 12 are constantly in tension, while the thick-walled lower part of the BHA 11 is partly in compression. The tension in the drill pipes avoids buckling of the drill pipe section. The torsional rigidity of the drill pipe section is, however, relatively small due to its slender construction.

(23) In practice, several types of speed controllers 20 have been developed and used, the control operation of which complies to a well-known PI controller, operable for providing a type of proportional action, P, and a type of integral action, I, or as PII controller, having a double integral action, for example. In the case of an electric drive system motor 18, for example, the speed controller 20 may be arranged to operate on a feedback from any or all of measuring variables such as the drive motor current, the rotational speed of the drive motor, and fluctuations in the drive motor current and rotational speed. This, for example, to control the energy flow in the drive system 15 by controlling any or both of these variables and to measure the rotational speed exerted by the drive motor 18 at the top end 14 of the drill string 12.

(24) The drive system 15 may operate in different modes, such as a so-called spinning mode and make-up mode, the present disclosure is primarily but not exclusively directed to the drill-mode, during which the driller aims to effectively grind or cut away material from an earth formation or geological formation by pushing and turning the drill bit 17 and flushing the borehole with drilling fluid or mud.

(25) The dynamic system comprised by the top drive, drill string and bottom hole assembly, can be represented by an or a combination of different model representations, among which, a state-space model, an equivalent electrical circuit model, an equivalent mechanical torsional spring-inertia model, a segmented model, a continuous-time model, a discrete-time model, a spectral or frequency domain model, and wave propagation models, such as mechanical and electrical type transmission line models.

(26) Referring to FIG. 2, for the purpose of explaining the present disclosure, there is shown a simplified computational dynamic model representation of the borehole drilling equipment 10 of FIG. 1 in the spectral domain. The spectral domain is represented by the complex argument or operator s, following the well-known Laplace transform theory.

(27) Reference numerals 34, 35, 36, 37, 38 and 39 in FIG. 2 denote summation operations, i.e. represented by small circles, having inputs, indicated by arrows pointing towards the circle, and an output, indicated by an arrow pointing away from the circle. The value at the output of a respective adder is the sum of the values at its inputs taken into to account the respective sign, i.e. plus or minus, indicated at a particular input in FIG. 2.

(28) The drill string 12 is modelled in accordance with a so-called ‘lumped string model’. In such a lumped string model, the drill string 12 as a whole, i.e. from its top end 14 to its bottom end 13, is simplified by a single torsional spring, experiencing a damping action caused by an internal and an external friction, operating along the drill string 12. Internal friction, among others, is caused by deformation losses inside the steel of the pipes of the drill string and is proportional to the wind-up/down speed of the drill string. Drilling fluid or drilling mud that is pumped down towards the BHA 11, for lubrication purposes, and cuttings from the drill bit 17 that are circulated back up the annulus, i.e. the space between the outer circumference of the drill string 12 and the borehole wall, to the earth surface, are also a cause of friction due to the speed difference of the drill string and the mud. All such drilling fluid or mud related friction forces are considered external friction.

(29) In the spectral domain model of FIG. 2, the drill string 12 is represented by its transfer function 40, and is assumed to have a torsion stiffness K.sub.s [Nm/rad], string inertia J.sub.s [kgm.sup.2], internal friction Y.sub.int [Nms/rad] and external friction Y.sub.ext [Nms/rad]. In the spectral domain, the torsion action 23 of the drill string 12 is modeled by its torsion stiffness K.sub.s and the Laplace transform operation 1/s, i.e. K.sub.s/s. The inertia action 24 is modeled by the reciprocal or inverse of the string inertia J.sub.s and the Laplace transform operation 1/s, i.e. 1/(sJ.sub.s). The internal and external friction actions Y.sub.int, Y.sub.ext operating along the drill string 12 are represented by elements 26, 27, respectively, together with their summation sign at a respective input of the summation operations 36, 37, respectively.

(30) The BHA 11 and drill bit 17 may be represented by a concentrated inertia J.sub.bha [kgm.sup.2] and concentrated friction Y.sub.bha [Nms/rad] at the bottom end 13 of the drill string 12. In the spectral domain, the concentrated inertia action 25 is represented by the reciprocal or inverse of J.sub.bha and the Laplace transform operation 1/s, i.e. 1/(sJ.sub.bha). The concentrated friction action Y.sub.bha is represented by element 28, together with its summation sign at a respective input of the summation operation 39.

(31) The rotational drive system 15 is represented by an inertia J.sub.td [kgm.sup.2] and its action is modelled by the reciprocal or inverse of the inertia J.sub.td and the Laplace transform operation 1/s, i.e. 1/(sJ.sub.td), represented by block 29.

(32) The rotational drive system 15 in combination with a variable frequency drive (VFD), if applicable, acting as a source of torque with high bandwidth is condensed in a transfer function D(s), represented by reference numeral 30. Reference numeral 31 represents a speed-torque transfer function A(s) of the speed controller 20. A speed filter 32 having transfer function G(s) may be internally present in or externally connect to an input of the speed controller 20, represented by summing operation 34, for inputting the measurement signal of the speed indicator or speed sensor 21, for example, measuring the top end drive speed ω.sub.td[rad/s]. Block 33 represents a torque filter having transfer function E(s) that may be present in the speed controller 20.

(33) In the model of FIG. 2, ω.sub.set [rad/s] represents the driller's reference speed, i.e. the rotational speed set by the driller at the input of the speed controller 20. The signal ω.sub.err represents the difference between the speed signal ω.sub.td representative of the actual rotational speed at the top end 14 of the drill string 12, and the reference or set rotational speed signal ω.sub.set and is formed by the summation operation 34, i.e. ω.sub.err=ω.sub.set−ω.sub.td, and is input to transfer function A(s) 31.

(34) The speed controller 20 operates, while drilling, to generate a commanded torque signal T.sub.ref [Nm] for controlling the drive system motor 18, represented by the transfer function D(s) 30. The actual torque provided by the drive system 15, i.e. at the output of the transfer function D(s) in the computational model of FIG. 2, is termed T.sub.act [Nm]. T.sub.meas represents the actual torque T.sub.act as indicated by the speed controller 20 through the filter transfer function E(s) of the torque filter 33, for example.

(35) The accumulated torque T.sub.acc [Nm] experienced by the drive system 15, i.e. block 29 in the model of FIG. 2, is represented by a difference between the torque T.sub.act applied by the drive system motor 18 and the actual load torque T.sub.s [Nm] at the top end 14 of the drill string 12, formed by the summation operation 35, i.e. T.sub.acc=T.sub.act−T.sub.s. The load torque operating at the drill bit 17 is represented by T.sub.bit [Nm] and the rotational speed of the drill bit, i.e. the BHA 11, is expressed by ω.sub.bit [rad/s].

(36) It is noted that in the system model shown in FIG. 2, the drill string 12 is represented by a two-port transfer function 40 in the spectral domain, the input variables of which are the bit torque T.sub.bit and the top end drive speed ω.sub.td, and the output variables of which are the load torque T.sub.s and the rotational bit speed ω.sub.bit. It is this choice of input and output parameters that explain the choice of the signs at the respective inputs of the summation operations 35, 37 and 38, as indicated.

(37) In the present description and claims a generalized two-port transfer matrix notation is used, in which the bit torque T.sub.bit and the bit speed ω.sub.bit at the downhole side form the dependent variables and the load torque T.sub.s and the top end drive speed ω.sub.td at the top end side form the independent variables.

(38) The lumped string model may be further simplified by replacing the string inertia J.sub.s by the concentrated inertia J.sub.bha and neglecting the concentrated friction Y.sub.bha.

(39) In use, at its top end 14, the drill string 12 can be pulled upwards by hoisting machinery, also called draw-works (not shown). On the bottom end 13, when on-bottom, the BHA 11 is resting at or touches the earth formation by the drill bit 17. This contrary to an off-bottom BHA 11 position, in which the drill bit 17 does not touch the earth formation at the bottom of the borehole. It will be appreciated that in off-bottom operation the damping or concentrated friction Y.sub.bha will be less than when the drill bit touches the earth formation at the bottom of the borehole, which translates into a much higher damping term than when off-bottom.

(40) The lumped string model does not take into account any travel time of torsional waves nor higher modes in the drill string. The transfer function 40 of the drill string 12 as a whole is, however, also representative of the transfer function of a fraction of the drill string. Accordingly, the drill string as whole may be represented in the spectral domain by cascading a plurality of two-port transfer functions, each representing a section or length of the drill string 12 from its top end 14 to its bottom end 13, as shown in FIG. 3.

(41) The distributed spectral domain string model 45 of FIG. 3 is comprised of N−1 two-port string transfer function sections or segments 45.sub.1, 45.sub.2, . . . , 45.sub.n, . . . 45.sub.N−1, with N being an integer 2. The dependent variables of section n, i.e. torque T.sub.n+1 and speed ω.sub.n+1 form the independent variables of the subsequent section n+1, i.e. torque T.sub.n and speed ω.sub.n, respectively, of this section n+1. The speed ω.sub.1 of the first section at the top end 14 of the drill string 12 is the top end drive speed ω.sub.td and the torque T.sub.1 of the first section is the actual load torque T.sub.s at the top end 14 of the drill string 12. The speed ω.sub.N of the last section N at the bottom end 13 of the drill string 12 represents the rotational bit speed ω.sub.bit while the torque T.sub.N operates at the BHA 11 and drill bit 17.

(42) The length x along the drill string 12, i.e. the distance in the hole to be drilled, is schematically indicated by vector 44, running from the top end 14 to the bottom end 13 of the drill string 12. In a particular section n of the drill string, the respective fraction or length of the drill string 12 is assumed to have a torsion stiffness K.sub.n,n+1, inertia J.sub.n,n+1, internal friction Y.sub.nint and external friction Y.sub.next.

(43) All elements are assumed to be linear, which is a valid assumption for identifying the effective ‘small-signal’ friction. Coulomb friction, i.e. a steady torque only depending on the direction of movement, or ‘Stribeck friction’, a torque that decreases with increasing speed, effectively acting as a destabilizing negative small-signal friction, which is regarded to be the driving force in persisting stick-slip oscillations during unfavorable settings of the top drive controller, are not represented with the linear model of FIG. 2. When Coulomb friction is to be taken into account, same may be estimated as a fixed value of, generally, a few Nm/meter or a few kNm/kilometer drill string length and may be added to the external mud related friction.

(44) As an alternative, the mechanically based computational model of the drill string of FIG. 3 may be represented by an electrical equivalent circuit diagram comprising electrical elements. FIG. 4 shows an electrical equivalent of a string fraction 45.sub.n of the distributed string model of FIG. 3.

(45) In the circuit diagram of FIG. 4 the respective fraction or length 45.sub.n of the drill string is modeled by an inductor L.sub.n with an inductance value L.sub.n=1/K.sub.n,n+1 [H]. The inertia of the respective fraction or length of the drill string is modeled by a capacitor C.sub.n with a capacitance value C.sub.n=J.sub.n,n+1 [F]. Internal and external friction are represented by a resistor R.sub.nint having a resistance value R.sub.nint [Ω] and connected in parallel to the inductor L.sub.n and a resistor R.sub.next having a resistance value R.sub.next [Ω] and connect in parallel to the capacitor C.sub.n, respectively. Note that admittance Y [A/V] and impedance or resistance R [V/A] are reciprocal, thus R.sub.nint=1/Y.sub.nint and R.sub.next=1/Y.sub.next.

(46) For the section n, current i.sub.n [A] and voltage u.sub.n [V] represent the independent variables torque T.sub.n and speed ω.sub.n, respectively, and current i.sub.n+1 [A] and voltage u.sub.n+1 [V] represent the dependent variables torque T.sub.n+1 and speed ω.sub.n+1, respectively.

(47) A realistic distributed string model may be built by using many dozens or even hundreds of fractions or lengths, each consisting of four parameters and two states (local speed and torque). Such a model resembles a ladder network used in transmission line models in electronics. The main difference, however, is that R.sub.nint being parallel to the inductor L.sub.n. In the electrical world the copper resistance of the electrical conductors is in series with the per length inductance, not in parallel as presented here in the mechanical string model of FIG. 4. Another difference is the often ignored resistance R.sub.next in electrical transmission line models as the electrical isolation between the conductors is mostly very high and its effect may be neglected although in coaxial cable transmission line models dielectric losses caused by this electrical isolation are considered. However, torsional viscous friction exerted by the mud and bore hole wall friction is very substantial in borehole drilling and have to be incorporated in the mechanical string model and, accordingly, in the electrical equivalent thereof. Hence, the standard electrical equivalent circuit diagram of a transmission line needs to be modified to be applicable to model a drill string, as shown by the ladder-network transmission line model of FIG. 5.

(48) The equivalent electrical circuit diagram of FIG. 5 is a cascade of the electrical equivalent circuit diagrams in accordance with FIG. 4 of the respective string fractions 45.sub.1, 45.sub.2, . . . , 45.sub.n, . . . 45.sub.N−1 of the distributed string model of FIG. 3 and may be represented by a two-port transfer function in the spectral domain. Input current i.sub.1 and input voltage u.sub.1, represented by voltage source 46, are the independent variables T.sub.1 and ω.sub.1, respectively, and output current i.sub.N+1, represented by current source 47, and output voltage u.sub.N+1 are the dependent variables T.sub.N and ω.sub.N respectively.

(49) A model such as shown in FIG. 5 would create a very large number of parameters to be identified. Fortunately the circuit in FIG. 5 may be written in the spectral or frequency domain for the limit case where N goes to infinity. It can be proven that for a drill string of length L and an infinite number N of string sections of length dx, such that N=L/dx.fwdarw.∞, a two-port transfer matrix in the spectral domain may be expressed as:

(50) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_1\\ T_1\end{array}\right)=\left[\begin{array}{cc}\cosh \left(\mathrm{yL}\right)& Z_0\sinh \left(\mathrm{yL}\right)\\ Y_0\sinh \left(\mathrm{yL}\right)& \cosh \left(\mathrm{yL}\right)\end{array}\right]\left(\begin{array}{c}{\omega }_N\\ T_N\end{array}\right)& \left(2\right)\end{array}$
with transfer matrix properties: γ=propagation coefficient [1/m] Y.sub.0=characteristic admittance [Nms/rad] Z.sub.0=1/Y.sub.0=characteristic impedance [rad/Nms] L=length of drill string [m] cos h(γL)=(e.sup.γL+e.sup.−γL)/2 sin h (γL)=(e.sup.γL−e.sup.−γL)/2 ω.sub.1=rotational speed at top end of the drill string [rad/s] T.sub.1=torque at the top end of the drill string [Nm] ω.sub.N=rotational speed at bottom end of the drill string [rad/s] T.sub.N=torque at the bottom end of the drill string [Nm]

(51) The propagation coefficient and characteristic admittance or impedance can be derived from the mechanical properties of the drill string. Assume a drill string constructed of steel pipes having specific mass ρ [kg/m.sup.3], shear modulus G [N/m.sup.2], specific viscous wall damping κ [Ns/m.sup.4], and specific viscous damping originating in material hysteresis δ [Ns/m.sup.2]. With s being the Laplace transform operator, in the spectral domain applies:
α=ρ.Math.s+κ [Ns/m.sup.4]  (3)
β=(G/s)+δ [Ns/m.sup.2]  (4)

(52) The dynamic properties of a length dx [m] of the drill string having a cross sectional polar moment IP=π/32.Math.(OD.sup.4−ID.sup.4) [m.sup.4], wherein OD=outer diameter [m] and ID=inner diameter [m] of a string pipe, can be expressed as:
α.Math.IP.Math.dx=ρ.Math.s.Math.IP+κ.Math.IP.Math.dx [Nms/rad]  (5)
β.Math.(IP/dx)=(G.Math.IP)/(s.Math.dx)+δ.Math.IP/dx [Nms/rad]  (6)
wherein: ρ.Math.s.Math.IP≡inertia J.sub.n,n+1, of the drill string section, κ.Math.IP.Math.dx≡wall damping or external friction Y.sub.next of the drill string section, (G.Math.IP)/(s.Math.dx)≡torsion stiffness K.sub.n,n+1 of the drill string section, and δ.Math.IP/dx≡material damping or internal friction Y.sub.nint of the drill string section.

(53) The resulting transfer matrix properties can be found from:

(54) $\begin{array}{cc}\mathrm{propagation}\mathrm{coefficient}y=\sqrt{\frac{\alpha .Math.\mathrm{IP}.Math.\mathrm{dx}}{\beta .Math.\mathrm{IP}\text{/}\mathrm{dx}}}=\sqrt{\frac{\alpha }{\beta }}.Math.\mathrm{dx}& \left(7\right)\\ \mathrm{characteristic}\mathrm{admittance}{Y}_0=\mathrm{IP}\sqrt{\alpha .Math.\beta }& \left(8\right)\\ \mathrm{characteristic}\mathrm{impedance}{Z}_0=1\text{/}{Y}_0=1\text{/}\left(\mathrm{IP}\sqrt{\alpha .Math.\beta }\right)& \left(9\right)\end{array}$

(55) Equation (2) can easily calculate, in the spectral domain, the behaviour of a uniform drill string of any length.

(56) In case the mechanical properties of the drill string and/or the internal and external friction are not uniform along the length of the drill string, same by modelled by a composite transfer matrix comprised of a matrix multiplication of a cascade of two-port transfer matrices in the spectral domain, the transfer matrix components of each respective transfer matrix of this cascade than represents the mechanical properties and friction or damping of the respective part or length of the drill string that may be assumed uniform, in accordance with the above transfer matrix (2).

(57) For example, assume the drill string 12 in FIG. 1 is modelled, in accordance with the electrical transmission line model described above, by a first section (a) extending from the top end 14 and a second section (b) terminating at the bottom end 13 in the borehole. Let ω.sub.1 and T.sub.1 the rotational speed and torque operating at the drill string at the top end, ω.sub.2 and T.sub.2 the rotational speed and torque operating at the interface of the first (a) and second sections (b) of the drill string, and ω.sub.3 and T.sub.3 the rotational speed and torque operating at the end of the drill string. With m.sub.a11, m.sub.a2, m.sub.a2, and m.sub.a22 representing the transfer matrix components of the transfer matrix of section (a) of the drill string and m.sub.b11, m.sub.b12, m.sub.b21, and m.sub.b22 representing the transfer matrix components of the transfer matrix of section (b) of the drill string, in accordance with equation (2) above.

(58) The two-port transfer matrix in the spectral domain for section (a) of the drill string may be expressed as:

(59) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_1\\ T_1\end{array}\right)=\left[\begin{array}{cc}{m}_{a11}& m_{a12}\\ m_{a21}& m_{a22}\end{array}\right]\left(\begin{array}{c}{\omega }_2\\ T_2\end{array}\right)& \left(10\right)\end{array}$

(60) The two-port transfer matrix in the spectral domain for section (b) of the drill string may be expressed as:

(61) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_2\\ T_2\end{array}\right)=\left[\begin{array}{cc}{m}_{b11}& m_{b12}\\ m_{b21}& m_{b22}\end{array}\right]\left(\begin{array}{c}{\omega }_3\\ T_3\end{array}\right)& \left(11\right)\end{array}$

(62) The two-port transfer matrix in the spectral domain of the drill string as whole is then found by substituting equation (11) in equation (10) resulting in:

(63) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_1\\ T_1\end{array}\right)=\left[\begin{array}{cc}{m}_{a11}& m_{a12}\\ m_{a21}& m_{a22}\end{array}\right].Math.\left[\begin{array}{cc}{m}_{b11}& m_{b12}\\ m_{b21}& m_{b22}\end{array}\right]\left(\begin{array}{c}{\omega }_3\\ T_3\end{array}\right)& \left(12\right)\end{array}$

(64) Performing the matrix multiplication provides:

(65) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_1\\ T_1\end{array}\right)=\left[\begin{array}{cc}{m}_{a11}.Math.m_{b11}+m_{a12}.Math.m_{b21}& m_{a11}.Math.m_{b12}+m_{a12}.Math.m_{b22}\\ m_{a21}.Math.m_{b11}+m_{a22}.Math.m_{b21}& m_{a21}.Math.m_{b12}+m_{a22}.Math.m_{b22}\end{array}\right]\left(\begin{array}{c}{\omega }_3\\ T_3\end{array}\right)& \left(13\right)\end{array}$

(66) It will be appreciated that the transfer matrix for a drill string 12 composed of more than two sections having different mechanical properties and/or experiencing different frictions, may be calculated likewise.

(67) A complete model of the borehole drilling equipment 10, i.e. speed controller 20, drive system 15, drill string 12 and BHA 11 and drill bit 17 can be derived by combining the two-port spectral domain transfer matrix representation (2) of the drill string 12 and respective two-port spectral domain transfer matrix representations of the BHA 11 and drill bit 17 on the one hand and of the speed controller 20 and drive system 15 on the other.

(68) For the BHA 11 and drill bit 17, the two-port spectral domain transfer matrix including friction experienced by the BHA 11 and drill bit 17 seen as a whole, i.e. concentrated, may be expressed by:

(69) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_N\\ T_N\end{array}\right)=\left[\begin{array}{cc}1& 0\\ s{J}_{\mathrm{bha}}+Y_{\mathrm{bha}}& 1\end{array}\right]\left(\begin{array}{c}{\omega }_{\mathrm{bit}}\\ T_{\mathrm{bit}}\end{array}\right)& \left(14\right)\end{array}$
wherein: J.sub.bha=concentrated bottom hole assembly inertia Y.sub.bha=concentrated bottom hole assembly friction ω.sub.bit=rotational bit speed ω.sub.N=rotational speed at bottom end of drill string T.sub.bit=torque operating at the drill bit T.sub.N=torque operating at bottom end of drill string s=Laplace transform operator.

(70) The two-port spectral domain transfer matrix for the rotational speed and torque operating at the top end 14 of the drill string 12 by the drive system 15 including the time delay τ.sub.d [s] caused by the speed controller 20 may be expressed by:

(71) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_{\mathrm{td}}\\ T_{\mathrm{ref}}\end{array}\right)=\left[\begin{array}{cc}1& 0\\ s{J}_{\mathrm{td}}.Math.e^{s{\tau }_d}& e^{s{\tau }_d}\end{array}\right]\left(\begin{array}{c}{\omega }_1\\ T_1\end{array}\right)& \left(15\right)\end{array}$
wherein: J.sub.t=top drive inertia τ.sub.d=speed controller time delay ω.sub.td=obtained rotational bit speed at top end of the drill string ω.sub.1=rotational speed at top end of the drill string T.sub.ref=commanded top end drive torque T.sub.1=torque at top end of drill string s=Laplace transform operator.

(72) Substituting equations (2) and (14) into equation (15) provides the two-port spectral domain transfer matrix Mbh of the borehole drilling equipment 10 as a whole, i.e.:

(73) 0$\begin{array}{cc}\left(\begin{array}{c}{\omega }_{\mathrm{td}}\\ T_{\mathrm{ref}}.Math.e^{-s{\tau }_d}\end{array}\right)=\left[\begin{array}{cc}{\mathrm{Mbh}}_{11}\left(s\right)& {\mathrm{Mbh}}_{12}\left(s\right)\\ {\mathrm{Mbh}}_{21}\left(s\right)& {\mathrm{Mbh}}_{22}\left(s\right)\end{array}\right]\left(\begin{array}{c}{\omega }_{\mathrm{bit}}\\ T_{\mathrm{bit}}\end{array}\right)& \left(16\right)\end{array}$
wherein Mbh.sub.11(s), Mbh.sub.12(s), Mbh.sub.21(s) and Mbh.sub.22(s) are the spectral transfer matrix components of the transfer matrix Mbh of the complete borehole drilling equipment 10 as indicated above, with:

(74) $\begin{array}{cc}\begin{array}{c}{\mathrm{Mbh}}_{11}\left(s\right)=\cosh \left(\mathrm{yL}\right)+Z_0\left(s{J}_{\mathrm{bha}}+Y_{\mathrm{bha}}\right)\sinh \left(\mathrm{yL}\right)\\ {\mathrm{Mbh}}_{12}\left(s\right)=Z_0\sinh \left(\mathrm{yL}\right)\\ {\mathrm{Mbh}}_{21}\left(s\right)=s{J}_{\mathrm{td}}\left\{\cosh \left(\mathrm{yL}\right)+Z_0\left(s{J}_{\mathrm{bha}}+Y_{\mathrm{bha}}\right)\sinh \left(\mathrm{yL}\right)\right\}-\\ \left\{Y_0\sin \left(\mathrm{yL}\right)+\left(s{J}_{\mathrm{bha}}+Y_{\mathrm{bha}}\right)\cosh \left(\mathrm{yL}\right)\right\}\\ {\mathrm{Mbh}}_{22}\left(s\right)=s{J}_{\mathrm{td}}{Z}_0\sinh \left(\mathrm{yL}\right)+\cosh \left(\mathrm{yL}\right).\end{array}\}& \left(17\right)\end{array}$

(75) Applying the rules of matrix inversion to equation (16) results in:

(76) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_{\mathrm{bit}}\\ T_{\mathrm{bit}}\end{array}\right)=\left[\begin{array}{cc}{\mathrm{Mbh}}_{22}\left(s\right)& -{\mathrm{Mbh}}_{12}\left(s\right)\\ -{\mathrm{Mbh}}_{21}\left(s\right)& {\mathrm{Mbh}}_{11}\left(s\right)\end{array}\right]\left(\begin{array}{c}{\omega }_{\mathrm{td}}\\ T_{\mathrm{ref}}.Math.e^{-s{\tau }_d}\end{array}\right)& \left(18\right)\end{array}$
since the determinant (Det) of the transfer matrix Mbh of the complete borehole drilling equipment equals 1, i.e. Det (Mbh)=Mbh.sub.11(s).Math.Mbh.sub.22(s)−Mbh.sub.12(s).Math.Mbh.sub.21(s)=1. In the notation of equation (18) the bit torque T.sub.bit and the bit speed ω.sub.bit at the downhole side form the dependent variables and the load commanded top end drive torque T.sub.ref and the top end drive speed ω.sub.td at the top end side form the independent variables.

(77) Observe that the rotational speed ω.sub.bit at the bottom end of the drill string and the torque T.sub.bit operating at the drill bit can be calculated, in the spectral domain, from equation (18) by applying matrix multiplication, such that:
ω.sub.bit=Mbh.sub.22(s).Math.ω.sub.td−Mbh.sub.12(s).Math.T.sub.ref.Math.e.sup.−sτ.sup.d  (19)
T.sub.bit=−Mbh.sub.21(s).Math.ω.sub.td+Mbh.sub.11(s).Math.T.sub.ref.Math.e.sup.−sτ.sup.d  (20)

(78) To ease and/or speed up the calculations, one may simplify equations (19) and (20) by, for example, assuming the concentrated friction Y.sub.bha=0 and Mbh.sub.11(s)=Mbh.sub.22(s).

(79) The spectral transfer matrix components of spectral the transfer matrix Mbh of the complete borehole drilling equipment 10 as indicated above are each frequency dependent transfer functions. To simplify the further description, the spectral transfer matrix components are written as Mbh.sub.11, Mbh.sub.12, Mbh.sub.21 and Mbh.sub.22.

(80) It has been observed that in practical systems, above a frequency threshold, due to the approximations and assumptions applied in accordance with equations (2)-(9), the requirement for matrix inversion, Det (Mbh)=Mbh.sub.11.Math.Mbh.sub.22−Mbh.sub.12−Mbh.sub.21=1 is no longer met. This has the effect of rising phases of the transfer matrix components Mbh.sub.11, Mbh.sub.12, Mbh.sub.21, and Mbh.sub.22 above this frequency threshold. Further, for higher frequencies, gains of the transfer matrix components rapidly increase due to the non-zero real part of the propagation coefficient γ in the functions cos h(γL)=(e.sup.γL+e.sup.−γL)/2 and sin h (γL)=(e.sup.γL−e.sup.−γL)/2.

(81) To compensate for these phenomena, the transfer matrix components Mbh.sub.11, Mbh.sub.12, Mbh.sub.21 and Mbh.sub.22 are modified by reducing the spectral bandwidth thereof such to include the fundamental mode and a limited number of higher order spectral modes of the transfer matrix components, and by adding, in the spectral domain, a time delay term to each of the transfer matrix components, the time delay being selected to enable a causal time domain solution.

(82) Suppose that the transfer function of a bandwidth reduction filter, i.e. a low-pass band filter, in the spectral domain is F(s) and that a time delay TF [s] is added to enable a causal time domain solution. Than, in the spectral domain, from equations (18) and (19), the rotational speed ω.sub.bit at the bottom end of the drill string and the torque T.sub.bit operating at the drill bit, applying the bandwidth reduction, i.e. the limitation or reduction of higher order spectral modes, and time delay operations, can be calculated from:
ω.sub.bit={Mbh.sub.22.Math.ω.sub.td−Mbh.sub.12.Math.T.sub.ref.Math.e.sup.−sτ.sup.d}.Math.F(s).Math.e.sup.−sτ.sup.F  (21)
T.sub.bit={−Mbh.sub.21.Math.ω.sub.td+Mbh.sub.11.Math.T.sub.ref.Math.e.sup.−sτ.sup.d}.Math.F(s).Math.e.sup.−sτ.sup.F  (22)

(83) In the time domain ω.sub.bit and T.sub.bit may be calculated from the inverse Laplace transform of equations (21) and (22), respectively.

(84) The bandwidth reduction may be obtained by any suitable filter function and can be based on, for example, at least one of the damping properties and length of the drill string, such to include in the calculation only higher order spectral modes having a spectral power density above a selected value and/or a maximum number of higher order spectral modes, such as 10 or preferably up to 20 or more preferably up to 30 spectral modes above the fundamental mode, for example.

(85) As low-pass band filter function F(s) an n.sup.th order low-pass Bessel filter operation may be used, wherein n≥4, preferably n≥10, more preferably n≥20, and tapering. For the purpose of alleviating the calculation of time-domain equivalents, a gradually decreasing tapering or window function is preferred, such as but not limited to, a cosine, a confined Gaussian, a Welch and a Lanczos window function or similar.

(86) Applying the filter function and delay multiplications in equations (21) and (22) above, with F(s) written as F, results in:
ω.sub.bit=(F.Math.e.sup.−sτ.sup.F.Math.Mbh.sub.22).Math.ω.sub.td−(F.Math.e.sup.−sτ.sup.F.Math.Mbh.sub.12).Math.T.sub.ref.Math.e.sup.−sτ.sup.d  (23)
T.sub.bit=−(F.Math.e.sup.−sτ.sup.F.Math.Mbh.sub.21).Math.ω.sub.td+(F.Math.e.sup.−sτ.sup.F.Math.Mbh.sub.11).Math.T.sub.ref.Math.e.sup.−sτ.sup.d  (24)
and written in matrix notification in accordance with equation (18):

(87) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_{\mathrm{bit}}\\ T_{\mathrm{bit}}\end{array}\right)=\left[\begin{array}{cc}{\mathrm{Mb}}_{11}& {\mathrm{Mb}}_{12}\\ {\mathrm{Mb}}_{21}& {\mathrm{Mb}}_{22}\end{array}\right]\left(\begin{array}{c}{\omega }_{\mathrm{td}}\\ T_{\mathrm{ref}}.Math.e^{-s{\tau }_d}\end{array}\right)& \left(25\right)\end{array}$
wherein Mb.sub.11, Mb.sub.12, Mb.sub.21 and Mb.sub.22 represent bandwidth limited and time delay term modified spectral transfer matrix components for the complete borehole drilling equipment including drill string inertia, borehole wall damping, stiffness and material damping, top end drive inertia, and bottom hole assembly inertia, in accordance with:

(88) $\begin{array}{cc}\begin{array}{c}{\mathrm{Mb}}_{11}=F.Math.e^{-s{\tau }_F}.Math.{\mathrm{Mbh}}_{22}\\ {\mathrm{Mb}}_{12}=-F.Math.e^{-s{\tau }_F}.Math.{\mathrm{Mbh}}_{12}\\ {\mathrm{Mb}}_{21}=-F.Math.e^{-s{\tau }_F}.Math.{\mathrm{Mbh}}_{21}\\ {\mathrm{Mb}}_{22}=F.Math.e^{-s{\tau }_F}.Math.{\mathrm{Mbh}}_{11}.\end{array}\}& \left(26\right)\end{array}$

(89) Equation (25) may be written in a generalized two-port spectral domain transfer matrix computational model for calculating down hole speed and/or down hole torque at any position x in the borehole or along the drill string, as schematically indicated by arrow 44 in FIG. 3, in accordance with:

(90) $\begin{array}{cc}\left(\begin{array}{c}{\omega }_{\mathrm{dh}}\\ T_{\mathrm{dh}}\end{array}\right)=\left[\begin{array}{cc}{M}_{11}\left(x\right)& M_{12}\left(x\right)\\ M_{21}\left(x\right)& M_{22}\left(x\right)\end{array}\right]\left(\begin{array}{c}{\omega }_{\mathrm{td}}\\ T_{\mathrm{ref}}.Math.e^{-s{\tau }_d}\end{array}\right)& \left(27\right)\end{array}$ wherein: M.sub.11(x), M.sub.12(x), M.sub.21(x) and M.sub.22(x)=bandwidth limited and time delay term modified spectral transfer matrix components at a particular position x along the length of the drill string, ω.sub.dh=down hole speed, T.sub.dh=down hole torque, ω.sub.td=obtained top end drive speed, T.sub.ref=commanded top end drive torque, τ.sub.d=speed controller time delay, and s=Laplace transform operator.

(91) Assuming a uniform mechanical composition of the drill string and uniform damping or friction properties along the drill string length, for a given position x [m] along the drill string or the borehole measured from a fixed reference position, i.e. the top end of the drill string, the respective modified transfer matrix components of equation (27) may be derived from equations (2-9), such that:

(92) $\begin{array}{cc}\begin{array}{c}{M}_{11}\left(x\right)=F.Math.e^{-s{\tau }_F}.Math.\left\{s{J}_{\mathrm{td}}{Z}_0\sinh \left(y.Math.x\right)+\cosh \left(y.Math.x\right)\right\}\\ M_{12}\left(x\right)-F.Math.e^{-s{\tau }_F}{Z}_0\sinh \left(y.Math.x\right)\\ M_{21}\left(x\right)=-F.Math.e^{-s{\tau }_F}.Math.\left\{s{J}_{\mathrm{td}}\cosh \left(y.Math.x\right)+Y_0\sinh \left(y.Math.x\right)\right\}\\ M_{22}\left(x\right)=F.Math.e^{-s{\tau }_F}.Math.\cosh \left(y.Math.x\right)\end{array}\}& \left(28\right)\end{array}$

(93) In the event of a non-uniform drill string and/or non-uniform or varying internal and external friction along the length of the drill string, same by modelled by a composite transfer matrix as disclosed above with reference to equations (10-13).

(94) Note that in equation (28) the modified spectral transfer matrix components are each spectral transfer functions and may be written in full notation as M.sub.11(x;s), M.sub.12(x;s), M.sub.21(x;s) and M.sub.22(x;s).

(95) In the time-domain, any or both of the down hole speed and down hole torque are efficiently calculated in a convolution operation processing the obtained rotational top end drive speed and top end drive torque and impulse response functions of the modified spectral transfer matrix components M.sub.11, M.sub.12, M.sub.21 and M.sub.22 over a period of time including one round trip time in the drill string. That is back and forth the drill string length.

(96) A convolution of two continuous functions g and h over a finite time range [−t.sub.0, t.sub.0] is generally expressed by:

(97) $\begin{array}{cc}\left[g*h\right]\left(t\right)=\underset{-t_0}{\overset{t_0}{\int }}g\left(u\right).Math.h\left(t-u\right)\mathrm{du}& \left(29\right)\end{array}$
and for discrete signals g′ and h′ over a number of samples n in the time domain over a finite time range [−p.sub.0, p.sub.0] by:

(98) $\begin{array}{cc}\left[g^{\prime }*h^{\prime }\right]\left(n\right)\equiv \underset{-p_0}{\overset{p_0}{.Math.}}{g}^{\prime }\left(v\right).Math.h^{\prime }\left(n-v\right)& \left(30\right)\end{array}$

(99) The impulse response functions of the modified spectral transfer matrix components M.sub.11, M.sub.12, M.sub.21 and M.sub.22 are calculated from inverse spectral to time domain transformation and can be suitably synthesized from commercially available computer software and only need to span a length of time at least equal to a single round trip time in the drill string which, in practice, equals a few seconds.

(100) Assume the impulse response functions of the modified spectral transfer matrix transfer function components M.sub.11(s), M.sub.12(s), M.sub.21(s) and M.sub.22(s) are expressed by H.sub.11(t), H.sub.12(t), H.sub.21(t) and H.sub.22(t), respectively. As the measured time varying top end drive speed ω.sub.td(t) and the commanded top end drive torque T.sub.ref(t) are only known from the past, the time-varying down hole speed ω.sub.dh(t) and/or down hole torque T.sub.dh(t) of the borehole drilling equipment are calculated by causal convolution operations, i.e. over a finite time range running from [−2t.sub.0, 0]:

(101) $\begin{array}{cc}{\omega }_{\mathrm{dh}}\left(t-t_0\right)=\underset{-2t_0}{\overset{0}{\int }}{\omega }_{\mathrm{td}}\left(u\right).Math.H_{11}\left(t-u\right)\mathrm{du}+\underset{-2t_0}{\overset{0}{\int }}{T}_{\mathrm{ref}}\left(u-{\tau }_d\right).Math.H_{12}\left(t-u\right)\mathrm{du}& \left(31\right)\\ T_{\mathrm{dh}}\left(t-t_0\right)=\underset{-2t_0}{\overset{0}{\int }}{\omega }_{\mathrm{td}}\left(u\right).Math.H_{21}\left(t-u\right)\mathrm{du}+\underset{-2t_0}{\overset{0}{\int }}{T}_{\mathrm{ref}}\left(u-{\tau }_d\right).Math.H_{22}\left(t-u\right)\mathrm{du}& \left(32\right)\end{array}$

(102) From equations (31) and (32) one can see that the total delay in the calculation amounts t.sub.0, i.e. the calculated down hole speed and down hole torque are delayed over t.sub.0, which is expressed as ω.sub.dh(t−t.sub.0) and T.sub.dh(t−t.sub.0), respectively. The convolution operations (31) and (32), which are calculated in practice based on discrete signals, are computationally very efficient in that the range [−2t.sub.0, 0] may equal a period of time including at least one round trip time in the drill string, in particular a period of time being a fraction longer, say 10-20% longer, than one round trip time. Thus, wherein to minimal amounts a one way travel time in the drill string.

(103) For calculating the bit speed and/or bit torque, synthesized impulse response functions of the modified spectral transfer matrix components Mb.sub.11, Mb.sub.12, Mb.sub.21 and Mb.sub.22 from equations (17) and (26) are to be applied, i.e. Hb.sub.11(t), Hb.sub.12(t), Hb.sub.21(t) and Hb.sub.22(t).

(104) FIGS. 6a, 6b, 6c and 6d show examples of a the impulse response functions Hb.sub.11(t), Hb.sub.12(t), Hb.sub.2(t) and Hb.sub.22(t) including drill string inertia, borehole wall damping, stiffness, material damping, top end drive inertia, bottom hole assembly inertia, bit friction and assuming a uniform drill strength and uniform damping, computed from a drilling equipment simulation test system developed by applicant. This simulation test system comprises an electric drive motor, acting as top drive and controlled by a speed controller, and a controllable electric load, mechanically coupled to the drive motor, simulating the drill string, BHA and drill bit.

(105) For verification purposes, the borehole drilling simulation test equipment has been operated while simulating a practical drilling operation with stick-slip and applying stick-slip mitigation control, showing the following properties:

(106) TABLE-US-00001 propagation coefficient γ = 3.1404 10.sup.−4 [1/m] characteristic admittance Y.sub.0 = 377.71 [Nms/rad] characteristic impedance Z.sub.0 = 0.0026476 [rad/Nms] top drive inertia J.sub.td = 1060 [kgm.sup.2] concentrated bottom hole inertia J.sub.bha = 105.5 [kgm.sup.2] concentrated friction Y.sub.bha = 4.03 [Nms/rad] drill string length L = 3804 [m] specific external friction κ = 398.2 [Ns/m.sup.4] specific internal friction δ = 2.757 10.sup.−8 [Ns/m.sup.2]

(107) Calculations in accordance with the present disclosure have been performed applying a tapered 25.sup.th order Bessel filter function F having a low-pass bandwidth of 5 Hz and with an added time delay term τ.sub.F of 2 seconds.

(108) In FIGS. 6a, 6b, 6c and 6d, time t−τ.sub.F [s] is plotted along the horizontal axis, i.e. the graphs 61-68 are depicted shifted over a time delay of 2 s.

(109) In FIGS. 6a, 6b, 6c and 6d, the solid graphs 61, 62, 63, 64 show the impulse response functions over a period of time t calculated for the drill string represented by a four segment or section model, in accordance with the model shown in FIG. 2. The dashed graphs 65, 66, 67, 68 represents the impulse functions over a period of time t calculated for a single segment or section drill string model, as illustrated in FIG. 3.

(110) FIG. 6a shows the impulse response function Hb.sub.11(t), 61, 65, which is dimensionless, representing the contribution of the top end drive speed ω.sub.td to the bit speed ω.sub.bit. FIG. 6b shows the impulse response function Hb.sub.12(t) [rad/Nms], 62, 66, representing the influence of the commanded top end drive torque T.sub.ref to the bit speed ω.sub.bit. FIG. 6c shows the impulse response function Hb.sub.21(t) [Nms/rad], 63, 67, illustrating the contribution of the top end drive speed ω.sub.td at the bit torque T.sub.bit. FIG. 6d shows the transfer of the commanded top end drive torque T.sub.ref to the bit torque T.sub.bit by the impulse response function Hb.sub.22(t), 64, 68, which is also dimensionless.

(111) From the impulse response functions shown, it can be derived that the top end drive speed ω.sub.td and the commanded top end drive torque T.sub.ref only contribute to the bit speed ω.sub.bit and bit torque T.sub.bit in the limited time period from about −1.5 s to 1.5 s and that the impulse response functions are zero or essentially zero outside this time period. More particular, within the time period of −1.5 s to 1.5 s the impulse response functions are also partly zero or essentially zero, i.e. from about −1 to 1 s. Accordingly, for a relatively long uniform drill string, convolution operations need to be performed over a very limited range only.

(112) The solid line 73 in FIG. 7c illustrates graphically, over a selected time period t, the calculated bit speed ω.sub.bit in revolutions per second [rpm], based on the four segment drill string impulse response functions in FIGS. 6a-6d and actual values of the top end drive speed ω.sub.td and the commanded top end drive torque T.sub.ref measured over time from the simulated drilling operation performed by the borehole drilling simulation equipment mentioned above in connection with FIGS. 6a-6d. The dashed line 74 illustrates the actual measured bit speed at the simulation test equipment.

(113) During the simulation, at time t=disclosure s form the start of the simulation, stick-slip mitigation control is applied, in accordance with the method and device disclosed in the International Patent Application WO2013/062409, filed in the name of applicant. In FIG. 7b, prior to t=disclosure s stick-slip occurs, and the speed of the drill bit fluctuates rapidly between zero and about 200 rpm. From the time of applying stick-slip mitigation control, the bit speed fluctuations rapidly decrease, as shown in FIG. 7b from time t=disclosure s.

(114) Graph 71 in FIG. 7b shows the calculated contribution to the bit speed ω.sub.bit by the top end drive speed ω.sub.td, and graph 72 in FIG. 7c shows how the commanded top end drive torque T.sub.ref contributes to the bit speed ω.sub.bit over a period of time t, both based on the four segment impulse response functions shown in FIG. 6a-6d. Graph 73 is the sum of both graphs 71 and 72.

(115) It is remarkably to note the accuracy of the calculated bit speed in accordance with the present disclosure and the measured values. FIG. 8 shows, to this end, an enlarged graphical presentation, illustrating the actual measured bit speed value, at the time interval directly on an enlarged tie scale, prior to and directly after applying stick-slip mitigation control at t=disclosure s. The dashed graph 74 represents the actual measured bit speed in the simulation and the solid line graph 73 represents the bit speed calculated in accordance with the present disclosure from a four segment drill string model. The dashed dotted graph 75 shows the bit speed calculated from a single segment model of the drill string, based on the drill string operation simulated as mentioned with reference to FIGS. 6a-6d.

(116) Both the four segment and single segment drill string computational model provide the same trend in the bit speed over time, in accordance with the measured values from the simulation test system.

(117) The contributions of the top end drive speed ω.sub.td and commanded top end drive torque T.sub.ref to the bit torque T.sub.bit calculated in accordance with the present disclosure for the measurement data obtained from the borehole drilling equipment simulation mentioned above with reference to FIGS. 6 and 7, are graphically illustrated in FIGS. 9a, 9b and 9c. Graph 76 in FIG. 9a illustrates the contribution of the top end drive speed ω.sub.td and graph 77 in FIG. 9b shows the contribution of the commanded torque T.sub.ref, both over a selected time period t, including the application of stick-slip mitigation control from t=disclosure s.

(118) The calculated bit torque shown in FIG. 9c is the sum of the contributions of FIGS. 9a and 9b. Graph 79 in FIG. 9d shows the measured bit torque T.sub.bit. For comparison reasons, the measured bit torque is displayed corrected for the time delay term τ.sub.F of 2 seconds used in the calculation of the bit torque. The difference in the torque values of graphs 78 and 79 results from Coulomb friction, i.e. a steady torque value that is not taken into account in the calculations performed. When Coulomb friction is to be taken into account, same may be estimated as a fixed value of, generally, a few Nm/meter or a few kNm/kilometer drill string length and may be added to the external mud related friction. In the present example, the Coulomb friction amounts 7.2 kNm. The accuracy and trend of the calculated value of the bit torque is remarkably.

(119) The graphs of FIGS. 7c and 9c show the bit speed and bit torque. However, as explained above, at each position down hole, i.e. along the drill string, down hole speed ω.sub.dh(t) and down hole torque T.sub.dh(t) may be calculated in accordance with the present disclosure, from synthesizing respective impulse response functions H.sub.11(t), H.sub.12(t), H.sub.21(t) and H.sub.22(t) of the modified spectral transfer matrix components M.sub.11(s), M.sub.12(s), M.sub.21(s) and M.sub.22(s).

(120) The positions are each time referred to a set fixed point of the borehole drilling equipment, such as the plane of the earth surface in which a borehole is drilled or another fixed reference position, such as the top end of the drill string. The impulse response functions for a particular position may be calculated in advance, i.e. in a first processing step, while the respective convolution operations are performed in a second processing step once the top end drive speed and top end drive torques are obtained.

(121) The down hole speed and down hole torque for the plurality of positions may be presented graphically, while drilling, at a display device in both a two-dimensional (2D) or three-dimensional or a spatial (3D) diagram as a function of time and position.

(122) FIGS. 10a and 10b illustrate in a 2D representation the estimated downhole speed as a function of the down hole position x, i.e. along the drill string, at a particular time t and four points in time earlier, each with 0.5 seconds difference, as indicated in the legend with the graphs. FIG. 10a represents the down hole speed at t=1811.5 s, that is prior to the application of stick-slip mitigation control at t=1840 s. FIG. 10b represents the down hole speed at t=1862 s, that is after the application of stick-slip mitigation control at t=1840 s. At the top of the graphs the momentary speed of the top drive is indicated, i.e. the speed at time t=1811.5 s and t=1840 s, respectively.

(123) In the 3D representation of FIG. 11, down hole speed is illustrated in function of the down hole distance x during a time interval t including the application of stick-slip mitigation control at t=1840 s, estimated, i.e. calculated in accordance with the present disclosure from top end measured data obtained of the simulation test system. FIG. 12 illustrates the down hole torque in function of the down hole distance x during a time interval t including the application of stick-slip mitigation control at t=1840 s, and calculated/estimated in accordance with the present disclosure from top end measured data obtained of the simulation test system

(124) It will be appreciated that such speed wave and/or down hole torque pattern, illustrated by FIGS. 10a, 10b, 11 and 12, may be displayed in near real time at a display of an operators console or the like, either on site or remote, such to continuously monitor the speed and torque acting at the drill string and the drill bit.

(125) It will be appreciated that the method according to the present disclosure is not only able to show and detect bit stick-slip, but is also able to detect stick-slip along the drill string, including higher mode stick-slip.

(126) As mentioned in the summary part, the relevant mechanical and string geometry operational parameters, for calculating the transfer matrix components may be obtained manually and/or digitally under operator control from a spread-sheet or table or the like. In a more sophisticated implementation, the method can be performed completely automatically and autonomously, while drilling commences, without the intervention of an operator or other skilled personnel, by acquiring the relevant mechanical and string geometry operational parameters from applying the teachings of the International Patent Application WO2014/098598, filed in the name of applicant, by applying at the commanded top en drive torque T.sub.ref a time-varying frequency torque control signal, such as illustrated in FIG. 2 by input signal T.sub.sweep, dashed arrow 42, and summation operation 41 at the input of the transfer function D(s) 30 and processing the obtained top end drive speed in an operational parameter control system 87 configured accordingly, see FIG. 14.

(127) FIG. 13 illustrates, by a simplified flow chart 50, the steps to be processed by a digital processor, computer or other data processing device, for estimating at least one of the down hole speed and/or down hole torque in accordance with an example of the present disclosure and presenting same, for example by displaying, or otherwise using such estimations for the purpose of drilling a borehole in an earth formation. The direction of flow is assumed from the top to the bottom of the sheet. Other directions are indicated by a respective arrow.

(128) As a first step, a suitable computational model representing the dynamic part including damping and friction of borehole drilling equipment 10 comprised by the drive system or top drive 15, drill string 12 and BHA 11 is selected, i.e. block 51 “Select two-port spectral domain transfer matrix computational model for borehole drilling equipment”, that is such that the down hole speed and down hole torque are the dependent variables and the rotational top end drive speed and top end drive torque are the independent variables. In general this computational model, such as the transmission line model disclosed above with reference to FIGS. 4 and 5, is selected once.

(129) Next, as shown by block 52, “Obtain transfer matrix components of the computational model selected”, for a particular distance downhole or a particular length along the drill string, the values of the transfer matrix components are obtained from the mechanical properties of the drill string 12, top drive 15 and BHA 11, including friction and damping values, for example as disclosed above with reference to the transmission line model and equations 2-9, 17 and 18.

(130) The transfer matrix components thus obtained are modified to admit a reduced number of higher spectral modes in the transfer matrix components, besides the fundamental mode. Further a time delay term is added to each transfer matrix component in the spectral domain, to enable causal time domain solutions. This as illustrated by block 55, “Modify the transfer matrix components by reducing spectral bandwidth to include the fundamental and a limited number of higher order spectral modes and add a time delay term to each of the spectral transfer matrix components”. Reference is made to the spectral transfer matrix components M.sub.11(s), M.sub.12(s), M.sub.21(s) and M.sub.22(s) and Mb.sub.11(s), Mb.sub.12(s), Mb.sub.21(s) and Mb.sub.22(s) detailed above.

(131) During drilling, with a particular sampling frequency to be selected, such as for example in the range of 4-20 ms, top end drive speed is obtained, such as measured by a speed sensor 21, concurrently with a value representative for the top end drive torque, i.e. block 54 “Obtain while drilling top end drive speed and top end drive torque”. However, any other suitable speed measurement technique or equipment may be used, such as sensorless estimating the top end drive speed from the voltage and current fed to an electric drive system motor 18. As explained above, the top end drive torque may be represented by the commanded top end drive torque set point T.sub.ref.

(132) Next, the down hole speed/bit speed and down hole torque/bit torque are calculated in the time domain, as illustrated by block 55, “Calculate down hole speed and/or down hole torque in the time domain from causal convolution operations of the obtained top end drive speed and top end drive torque and impulse responses of the modified transfer matrix components”. The impulse responses are obtained by inverse spectral to time domain transformation which may be performed, by the computer using an Inverse Fast Fourier Transformation, IFFT, algorithm, for example. Calculating the down hole speed/bit speed and/or down hole torque/bit torque from a convolution operation processing is advantageous in view of the short range over which the impulse response functions of the matrix components contribute.

(133) The thus computed speeds and torques may be presented in near real time to a driller or supervisor or other body, e.g. by displaying same, and/or recorded for further analysis, i.e. block 56, “Display calculated down hole speed and/or down hole torque for a plurality of down hole positions while drilling”.

(134) When in the course of a drilling operation the length of the drill string 12 is increased, for example, or if the drilling environment significantly changes, decision block 57, “Drilling environment/equipment changed?”, outcome “Yes”, the process as illustrated by blocks 52-56 is repeated to obtain adjusted transfer matrix components, if applicable. Otherwise, decision block 59 result “No”, the calculation and presentation of the down hole speed and/or down hole torques continues from the obtained top end drive speed and torque.

(135) It will be appreciated that the inertia J.sub.td of the drive system 15 and the inertia J.sub.bha of the bottom hole assembly 11 will not have to be established anew each time the drill string is extended. The inertia J.sub.td of the drive system may also be known beforehand from data provided by the operator of the borehole drilling equipment, for example, such to considerable ease the process. Adjusted transfer matrix components may be obtained as disclosed above.

(136) As mentioned in the summary part above, for computing the transfer matrix components of a selected computation model and their time-domain equivalents and impulse response functions, spectral domain transformations and inverse or spectral to time domain transformations, matrix manipulations, convolution operations and filter operations, the computer may be equipped with commercially available computer software such as MATLAB™ or any other commercially available computer program for computing spectral domain transformations and inverse transformations of transfer functions and impulse response functions of a dynamic system, matrix manipulations, convolution operations and filter operations, as well as proprietary discrete time modelled software programs written, for example, in C, C+ or C++.

(137) FIG. 14 schematically shows a device 80 for automatically estimating down hole speed and/or down hole torque of borehole drilling equipment 10, operated either on-bottom or off-bottom, while drilling a borehole. In addition to the borehole drilling equipment 10 shown in FIG. 1, a computer controlled drilling operation observer system 81 is provided. The observer system 81 comprises a computer or electronic processing device 82, an input interface 83, such as keyboard, a touch screen or the like, for selecting a representative two-port spectral domain transfer matrix computational model of the borehole drilling equipment 10 and/or for inputting operational parameter values, such as tuning parameters of the speed controller 20, if any, of the borehole drilling equipment 10 for determining transfer matrix components of the computational model. Speeds and torques calculated by the processing device 82 may be provided at an output interface 84, such as a graphical display, a printer or plotter, or a data evaluation module for evaluating the results obtained.

(138) The computational model, operational and tuning parameters, filter functions, time delay terms and other values in accordance with the present disclosure, to be used by the observer system 81, may be stored at and retrieved from a database 85, accessible from the observer system 81. The database 85 may be remote from the observer system 81 and connected by a communication network 86, for example, or may be comprised by a medium, including a transitory and non-transitory medium, that can be read by the computer or processing device 82.

(139) The rotational top end speed of the drill string 12, for example measured at its top end 14 by a speed indicator or speed sensor 21, is provided via input 89 to the processing device 82. The top end drive torque may be obtained from the speed controller 20 through a data input/output 88 for providing tuning parameters and/or feedback to the observer system 81.

(140) The processing device 82 is programmed, or may be programmed from a computer program product having a medium storing computer program code that can be read by the processing device 82, to calculate an estimated down hole speed and/or down hole torque from the obtained rotational top end drive speed and top end drive torque in accordance with the present disclosure. Further, the processing device 82 may be arranged to operate automatically and autonomously for obtaining the operational parameter values required for calculating the transfer matrix components by implementing an operational parameter control system in accordance with Patent Application WO2014/098598, filed in the name of applicant, which is herein incorporated by reference, and schematically indicated by reference numeral 87.

(141) The observer system 81 may be arranged or disposed remote from the borehole drilling equipment 10, such that the connections 88, 89 with the speed controller 20 may be arranged as a single data network connection or a connection with a telecommunication network, either wired or wireless.

(142) The speed controller 20 and/or a separate electronic controller and/or the system 81 may be further arranged for applying stick-slip mitigation control disclosed in the International Patent Application WO2013/062409, filed in the name of applicant, which is herein incorporated by reference.

(143) The present invention and disclosure is not limited to the embodiments as disclosed above, and can be modified and enhanced by those skilled in the art within the scope of the present disclosure as disclosed in the appended claims, without having to apply inventive skills.