Published 6 issues per year
ISSN Print: 2152-5080
ISSN Online: 2152-5099
IF:
1.7
5-Year IF:
1.9
Immediacy Index:
0.5
Eigenfactor:
0.0007
JCI:
0.5
SJR:
0.584
SNIP:
0.676
CiteScore™::
3
H-Index:
25
Indexed in
STOCHASTIC DRILL-STRING DYNAMICS WITH UNCERTAINTY ON THE IMPOSED SPEED AND ON THE BIT-ROCK PARAMETERS
Abstract
A drill string is a slender structure with nonlinear dynamics; it is an equipment used in the oil industry to drill rock in the search of oil and gas. The aim of this paper is to model the uncertainties related to the speed imposed at the top and uncertainties related to the bit-rock parameters, and to investigate how these uncertainties propagate throughout the system. The continuum system is linearized about the prestressed configuration, the finite-element model is applied to discretize the system, and then a reduced-order model is constructed using normal modes of the linearized system; only torsional and axial vibrations are considered in the analysis. A constant rotational speed is imposed at the top and a nonlinear bit-rock interaction acts at the bottom. A probabilistic approach is used to model the uncertainties and the Monte Carlo method is used to approximate the stochastic differential equations.
KEYWORDS: drill-string dynamics, stochastic dynamics, uncertainty quantification, bit-rock interaction
1. INTRODUCTION
The oil and gas industry has an expressive weight in the world economy. One important step involved in the exploitation of oil is the drilling process. In this sense, one should control the drill-string dynamics to avoid accidents and waste of resources. A general vibration perspective overview of the process of oil well rotary drilling can be found in the works of Spanos and co-workers [1, 2].
The drill-string dynamics is complex and has many sources of uncertainties. There are uncertainties related to the model as well as to the parameters of the model; here, model refers to the structure model, fluid-structure interaction model, bit-rock interaction model, and so on. Computational models must take into account these uncertainties to improve their predictions.
There are few articles concerned with uncertainty quantification and stochastic analysis of the drill-string dynamics; in particular, we may cite [3–8]. Spanos et al. [3] analyzed stochastic lateral forces at the bit. Kotsonis and Spanos [4] analyze a random weight-on-bit using a simple two degrees of freedom drill-string model. Ritto and co-workers proposed a probabilistic model for the bit-rock interaction model [5] using the nonparametric probabilistic approach [9, 10]. Then, they proposed a stochastic identification procedure to compute the dispersion parameter of the probabilistic model [6]. Later, they proposed a robust optimization problem to maximize the rate of penetration of the system without damaging it [7]. Finally, they analyzed a random weight-on-hook [8].
Some computational models have been developed to analyze the coupling between two or three vibration directions (see, for instance, the works of Yigit and Christoforou [11–13], or Khulief et al. [14], and also Sampaio and co-workers [15, 16]). A more complete model with coupled axial, torsional, and lateral vibrations has been developed by Tucker and Wang [17], where the Cosserat theory is applied, and by Ritto et al. [5], where the nonlinear Timoshenko beam theory is applied.
The aim of this paper is to analyze the stochastic drill-string dynamics with uncertainty on the top speed and on the bit-rock parameters, and to investigate how these uncertainties propagate throughout the nonlinear dynamical model. In this paper, coupled axial and torsional displacements are considered in the analysis, and the bit-rock interaction model used is based on [18], as done in [13, 14, 16], for instance.
The paper is organized as follows. In Section 2 the deterministic model is depicted in detail. In Section 3 the probabilistic model is constructed for the uncertain parameters, using the maximum entropy principle. Finally, in Section 4 the numerical results are analyzed and in Section 5 concluding remarks are made.
2. DETERMINISTIC SYSTEM
In the present model only the axial and the torsional vibrations of the column are considered. In order to focus attention on the uncertainty analysis (which is the goal of this paper), the lateral vibrations are assumed to be sufficiently small. The equation of motion for the axial vibration is given by
|
(1) |
where u is the axial displacement, ρ is the density of the column, A is the cross-sectional area, E is the elasticity modulus, ƒ is the force per unit length, and ƒNL is related to nonlinear and coupling terms. The equation of motion for the torsional vibration is given by
|
(2) |
where θx is the angular rotation about the x-axis, Ip is the cross sectional polar moment of inertia, G is the shear modulus, tor is the torque per unit length, and torNL is related to nonlinear and coupling terms.
In our modeling we assume finite strains, which consequently makes the equations of motion coupled and nonlinear. Therefore, it is easier to derive the nonlinear terms of the equations of motion applying the extended Hamilton principle. Defining the functional Π by
|
(3) |
where U is the potential strain energy, T is the kinetic energy, and W is the work done by the nonconservative forces and any force not accounted for in the potential energy. The first variation of Π must vanish:
|
(4) |
2.1 Some Definitions
Some definitions are necessary for the subsequent sections. Let X be a position in the initial configuration and x a position in the deformed configuration. The displacement field d = (ux , uy , uz ) written on the nondeformed configuration is given by (x − X):
|
(5) |
The velocity field v is the material time derivative of the displacement field and is written as
|
(6) |
where the material time derivative (d/dt) is denoted by a superposed dot. The stress-strain relationship (written in Voigt notation) is given by S = [D]ϵ, where S is the second Piola–Kirchhoff stress tensor S = (σx , τxy , τxz ), ϵ is the Green-Lagrange strain tensor ϵ = (ϵx , 2ϵxy , 2ϵxz ), and
|
(7) |
If finite strains are considered:
|
(8) |
Note that the strain tensor [E] might be derived using [E] = 1/2([F]T [F] − [I]), where [F] = dx/dX is the deformation gradient tensor and [I] is the identity operator, or using [E] = 1/2[dd/dX + (dd/dX)T + (dd/dX)T dd/dX].
2.2 Finite-Element Discretization
The finite-element model is constructed using two-node elements with two degrees of freedom per node (axial displacement u and angular rotation θx). The finite-element approximation of the displacement fields are then written as u(ξ, t) = Nu(ξ)ue(t) and θx(ξ, t) = Nθx(ξ)ue(t), where ξ = x/le is the element coordinate, N are linear shape functions, and element displacement is
|
(9) |
where exponent T means transposition.
2.3 Kinetic Energy
The kinetic energy is written as
|
(10) |
where V is the region of integration and L is the length of the column. The polar moment of inertia is written as Ip = ∫A(y2 + z2)dA, and the above expression was simplified by the approximation cos θx ≈ 1 and sin θx ≈ 0. The first variation of the kinetic energy, after integrating by parts in time, may be written as
|
(11) |
which yields the constant mass matrix [M]. The element mass matrix is written as:
|
(12) |
2.4 Strain Energy
The strain energy is given by
|
(13) |
Substituting the constitutive equation S = [D]ϵ and computing the first variation of the strain energy yield
|
(14) |
Following the analysis done in [7], the element stiffness matrix is written as
|
(15) |
where the space derivative (d/dξ) is denoted by (′). The element geometric stiffness matrix is written as
|
(16) |
where u′ = N′u ue / le , θ′x = N′θu ue / le , I22 = ∫A(y2 + z2)dA and Ip4 = ∫A(y4 + z4)dA.
2.5 Gravity
The work done by gravity is written as
|
(17) |
where g is the gravity acceleration. The variation of Eq. (17) gives
|
(18) |
and the discretization by means of the finite-element method yields the force element vector
|
(19) |
2.6 Bit-Rock Interaction
The bit-rock interaction couples the axial and torsional vibrations [13, 14, 16]; the bit-rock forces are semi-empirical expressions that are added a posteriori to the model. The torque applied at the bit is given by:
|
(20) |
where ωbit is the rotational speed of the bit, μbit is a factor that depends on the bit cutting characteristics, α1 and α2 are constants that depend on the rock properties, and ƒbit(t) = ƒc + ƒa(t) is the weight-on-bit. The weight-on-bit is the sum of the initial (constant) reaction force at the bit ƒc and a variable force that is time-dependent ƒa(t). This variable force is modeled by harmonic force [19], ƒa(t) = ƒ0 sin(ω0t) + 0.01 ƒ0, where the constant term on the right is added to avoid the possibility of bit bounce (not considered in the model). Note that this axial force does not depend on the response of the system, which makes the axial vibration much simpler than the torsional one.
2.7 Initial Prestressed Configuration
Before starting the rotation about the x axis, the column is put down through the channel until it reaches the soil. At this point, the forces acting on the structure are: the initial reaction force at the bit, the weight of the column, and the supporting force at the top. In this equilibrium configuration, the column is prestressed. Above the neutral point the structure is tensioned and below is compressed. If the reaction force increases, the neutral point moves up, increasing the length of the compressed part. To calculate the initial prestressed state, the column is clamped at the top and, consequently,
|
(21) |
where fg is the force induced by the gravity and fc is the vector related to the initial reaction force at the bit ƒc. (Note that ƒc = [0 0 … ƒc 0]T.)
2.8 Final Computational Dynamical Model
After computing the prestressed configuration (Section 2.7), the resulting displacement is used to obtain [Kg(uS)] (constant matrix). Finally, we assume small motions about this configuration; i.e., [Kg(u)] ≈ [Kg(uS)]. Introducing u = u − uS, the computational dynamical model can then be written as
|
(22) |
in which [M] and [K] are the mass and stiffness matrices. The proportional damping matrix [C] = α[M] + β([K] + [Kg(uS)]) (α and β are positive constants) is added a posteriori to the computational model. The initial conditions are defined by u0 and v0. The force vector related to the bit-rock interaction (force and torque at the bit) is fbit and the imposed rotation at the top (Dirichlet boundary condition) is included in g, which is the source force (input energy).
2.9 Reduced Computational Model
To speed up the numerical simulations, the computational model is reduced projecting the nonlinear dynamical equation on a subspace spanned by an appropriated basis. In the present paper, the basis used is made up of suitable normal modes. The normal modes are constructed from the following generalized eigenvalue problem:
|
(23) |
where φi is the ith normal mode and ωi is the corresponding natural frequency. The reduced model is written as
|
(24) |
in which q0 and v0 are the initial conditions and where [Φ] is the (m × n) real matrix composed of the n selected normal modes and
|
(25) |
are the reduced matrices.
3. PROBABILISTIC MODEL
The maximum entropy principle [20] is used to construct the probabilistic models for the uncertain parameters. This is, in the opinion of the authors, the best strategy to choose a probabilistic model (when data are limited). The resulting distribution obtained by the maximum entropy principle is always consistent with the information/constraints given about the random variable. If the information is consistent with the physics, then so are the obtained probability density functions.
The parameters modeled as random variables are the parameters of the bit-rock interaction model [ƒc, μbit, α1, α2; see Eq. (20)] and the imposed rotational speed at the top ω0. The reason to choose these parameters is that the parameters of the bit-rock interaction are hard to identify, and the top speed is the main driven force of the whole system; therefore, it is also a quantity of interest. In a real drilling operation, this speed might be different from its nominal value; it may vary, or it may be a little lower or upper, compared to its nominal value.
Let {Y1, Y2, Y3, Y4, Y5} represent the random variables related to the parameters {ƒc, μbit, α1, ω0}. The probability density functions of the random variables are obtained solving the following optimization problem:
|
(26) |
where p∗Y is the optimal probability density function such that ∀pY ∈ 𝒞p, S(p∗Y) ≥ S(pY), and 𝒞p is the set of admissible probability density functions that respect the available information presented in the sequence. The Shannon entropy measure S is defined as [21]:
|
(27) |
where supp is the support of the probability distribution. It should be noticed that if no information is known about the dependency among the random variables, the Maximum Entropy Principle yields independent random variables; which is the case of the present analysis. We will reason the available information for all of the random variables (Y) the following way:
- To keep the consistency with the formulation used, all the random variables should be positive (ƒc is a particular case, see the note below). Therefore, we consider the support ]0, +∞[.
- As an initial computational model, we will take the nominal parameters to be the mean of the random variables. Therefore, the mean value is given ℰ{Y} = Y.
- We would like the random variables to have small probabilities as Y approaches zero, since the probability of Y to be close to zero should be very small. To apply this condition we use ℰ{ln(Y)} = |c| < +∞.
Note: the random variable Y1 (related to ƒc) is always positive, but the reaction force should be negative. As a trick, we consider (−Y1) in the simulations; hence, the same reasoning applies. Note that we do not allow negative values for the rotation speed (ω0) and for the bit-rock parameters (μbit, α1, α2) to be consistent with our formulation. Nevertheless, this was a choice, the support could be set more flexibly or differently, depending on how one reasons things. Thus, there is more than one way to reason the available information; for instance, if we only give the information of the supports [Ymin, Ymax] (bounded support) this yields the uniform probability density function in this range.
Applying the maximum entropy principle with the information given in the above list (and, of course, the normalization condition) yields the gamma probability density function with mean Y and coefficient of variation δY = σY / Y, where σY is the standard deviation of Y. Finally, verify that this random variable has finite variance, which makes sense from the physical point of view. In addition, this random variable verifies that the random torque at the bit Tbit(t) and the random rotational speed of the bit Ωbit(t) are second-order random processes; i.e., they have finite variance, which makes sense from the physical point of view. As a final comment, this probabilistic model is not necessarily the best one; it is simply the least prejudice probabilistic model given the above constraints. Yet, we emphasize that the procedure of constructing a probabilistic model using the maximum entropy principle is preferable than just picking up an ad hoc probabilistic model; on the other hand, if data are available, there are other interesting strategies, such as adjusting the coefficient of a polynomial chaos expansion, Bayesian statistics, etc.
Finally, the deterministic system [Eq. (28)] becomes stochastic:
|
(28) |
where U(t) is the random response. The source term G(t) is random because of Y5 (related to the imposed speed at the top ω0), and the bit-rock interaction force vector Fbit is random because of Y1, Y2, Y3, and Y4 (related to the parameters of the bit-rock interaction model ƒc, μbit, α1, and α2). Besides these two random vectors, the random character of Y1 (related to the initial reaction force at the bit ƒc) makes the initial prestressed configuration [Eq. (21)] and the natural frequencies and normal modes [Eq. (23)] random. Therefore, [ψ], [Mr], [Cr], and [Kr] are random matrices.
The Monte Carlo method [22] is used to obtain the statistics of the responses of Eq. (28); i.e., for each realization of the independent random variables Yi (i = {1, 2, … 5}) the deterministic system is numerically integrated. The convergence is computed by the function conv(ns) = (1 / ns) ∑nsj=1 ∫ t1t0 ‖U (t, sj)‖2 dt, where ns is the number of Monte Carlo samplings.
4. NUMERICAL RESULTS
First the deterministic system is simulated and, then, each Yi is considered separately for δYi = {5%, 20%} (i = {1, 2, …, 5}) (see Section 3). The influence of each random parameter in the response of the system (axial and torsional displacements) is analyzed. Finally, all the random parameters are taken into account together. For the numerical integration in time, the explicit Runge-Kutta scheme of fourth-order and adapted time step was applied.
Typical values of a drill string are used: the 1600 m drill-string is composed of steel tubes (E = 200 GPa, ρ = 7850 kg/m3). The first 1400 m are composed of tubes with external diameter 120 mm and internal diameter 9.5 mm; the last 200 m are composed of thicker tubes (bottom hole assembly) with external diameter 150 mm. The parameters of the bit-rock interaction model are ƒ0 = 50 kN, ƒc = −100 kN, α1 = 2, α2 = 1, and μbit = 0.04 (data used in [14]). The damping constants are α = 1 × 10−1 and β = 8 × 10−5, and we set [Cr(1, 1)] = 0, which means that there is no damping related to the axial rigid body mode. The column is discretized in 120 finite elements, and the reduced-order model is composed of 10 torsional and 2 axial modes (one is the rigid body mode).
4.1 Deterministic Response
In this section, we analyze the relevance of considering the prestressed state for this problem (Section 4.1.1), and we investigate how many axial and torsional modes are needed in the reduced-order model (Section 4.1.2).
4.1.1 Natural Frequencies
Tables 1 and 2 show a comparison of the torsional and axial natural frequencies of the system with and without the prestressed configuration. The percent errors of the natural frequencies are very small; a similar result is verified for the normal modes (small errors). The conclusion is that the prestressed configuration is not important for this problem (at least for the parameters used in the analysis). It happens that it affects considerably the lateral response, which is not taken into consideration in the present analysis. Despite this fact, the prestressed configuration is taken into account in the numerical analysis, since it does not increase the computational time.
TABLE 1: First torsional natural frequencies in Hertz; comparing the system with and without the prestressed configuration.
| Configuration | First torsional natural frequencies (Hz) | |||||||||
| No prestress | 0.3679 | 1.2275 | 2.2039 | 3.2179 | 4.2411 | 5.2575 | 6.2380 | 7.1004 | 7.8407 | 8.7130 |
| Prestressed | 0.3681 | 1.2281 | 2.2049 | 3.2194 | 4.2431 | 5.2599 | 6.2407 | 7.1028 | 7.8435 | 8.7168 |
| Percent error | 0.06% | 0.05% | 0.05% | 0.05% | 0.05% | 0.05% | 0.04% | 0.03% | 0.04% | 0.04% |
TABLE 2: First axial natural frequencies in Hertz; comparing the system with and without the prestressed configuration.
| Configuration | First axial natural frequencies (Hz) | |||||||||
| No prestress | 1.4168 | 3.0166 | 4.7244 | 6.4701 | 8.2178 | 9.9314 | 11.5382 | 12.9613 | 14.3927 | 16.0154 |
| Prestressed | 1.4173 | 3.0179 | 4.7266 | 6.4731 | 8.2216 | 9.9358 | 11.5428 | 12.9659 | 14.3985 | 16.0226 |
| Percent error | 0.04% | 0.04% | 0.05% | 0.05% | 0.05% | 0.04% | 0.04% | 0.04% | 0.04% | 0.04% |
Figure 1 shows the frequency spectrum of the axial speed of the bit (rate of penetration). For the dynamics analyzed, only one axial mode (together with a rigid body mode) is necessary to get a good result; it should be noticed that in the axial direction the deterministic system is linear but the stochastic system is nonlinear due to the relation sin(Ω0), where Ω0 is the random variable related to the imposed top speed. Things are more complicated when the torsional vibration is analyzed because the bit-rock interaction is nonlinear. Figure 2 shows a comparison of the response of the rotational speed of the bit when one torsional mode is used and when 10 torsional modes are used in the reduced-order model. It is clear that one torsional mode is not sufficient to represent the dynamics of this system. Note also the stick-slip behavior of the response; i.e., sometimes the bit gets stuck (ωbit = 0), then it slips, with rotational speed varying from 0 to 300 rpm (the imposed nominal speed at the top is 100 rpm).
FIG. 1: Frequency spectrum of the rate of penetration of the bit for two approximations: 10 torsional modes and 10 axial modes versus 10 torsional modes and 1 axial mode. The two peaks observed are related to the first axial natural frequency (1.41 Hz) and to the imposed rotation at the top (100 rpm = 1.67 Hz).
FIG. 2: Angular speed of the bit for two approximations, 10 torsional modes and 10 axial modes versus 1 torsional mode and 10 axial modes: (a) time response and (b) zoom image.
4.2 Stochastic Response
A careful analysis was carried out separately for each random variable and then together for all random variables. Figure 3 shows the convergence function for the Monte Carlo samplings with δYi = 0.20 (i = 1, 2, …, 5); a reasonable convergence is achieved with 400 samplings.
FIG. 3: Convergence of the Monte Carlo simulation.
Analyzing the sensitivity of the bit angular speed for each case, a curious result is observed in which the random variables considered in the analysis affect the dispersion of the response (bit angular speed) in a similar manner: the dispersion of the response (standard deviation over the mean at each instant) is about 0.85 [see Fig. 4(a)], regardless of the uncertainty analyzed. This happens because in all Monte Carlo samplings there is stick-slip, and the bit angular speed oscillates between 0 and 300 rpm. However, the dispersion in the frequency domain is quite different [see Fig. 4(b)].
FIG. 4: Dispersion (standard deviation over the mean at each instant) ofthe bit angular speed for δY5 = 0.05 and δYi = 0.20 (i = 1, 2, …, 5): (a) time domain and (b) frequency domain.
It should be noted that the level of uncertainty of the random variables Yi (i = 1, 2, …, 5) was assumed to be δ = 0.05 or δ = 0.20. Actually, a stochastic inverse problem has to be performed to identify the dispersion level of each parameter of the system [for instance, we would expect the uncertainty related to the top speed (ω0) to be smaller than the uncertainty of a bit-rock interaction parameter such as α1]. The identification (which is a stochastic inverse problem) might be done, for example, using the maximum likelihood method [6] or using a Bayesian approach [23].
We are mostly concerned with the torsional vibration of the column at the bit (bit angular speed); nevertheless, the bit axial speed (which is the rate of penetration) is also analyzed. It should be mentioned that for the model considered in the present analysis the axial vibration of the column is not affected by uncertainties in the parameters μbit, α1 and α2, and is little affected by uncertainties in the parameter ƒc. Note that the axial vibration is coupled with the torsional one only through the rotational speed at the top (see Section 2.6), and the initial prestressed configuration (related to ƒc) does not have much influence on the axial vibration (see Section 4.1.1).
Figures 5–7 help us to analyze how the uncertainties are propagate throughout the nonlinear system. Figure 5 shows the 95% confidence region of the random response (axial and angular speed of the bit) together with the response of the nominal model, considering an uncertain top speed with δY5 = 0.05; and Fig. 6 shows the same graphic when all the parameters are considered uncertain: δYi = 0.20 (i = 1, 2, …, 5). As expected, as the uncertainties on the random input parameters increase, the uncertainty of the response also increases.
FIG. 5: (a) Axial speed and (b) angular speed of the bit: 95% confidence region (thick dark line), mean of the stochastic response (thin dark line), and response of the nominal model (blue line), considering an uncertain top speed with δY5 = 0.05.
FIG. 6: (a) Axial speed and (b) angular speed of the bit: 95% confidence region (thick dark line), mean of the stochastic response (thin dark line), and response of the nominal model (blue line), when all the parameters are considered uncertain with δYi = 0.20 (i = 1, 2, …, 5).
FIG. 7: Frequency spectrum of the bit angular speed (a) δY5 = 0.05 and (b) δYi = 0.20 (i = 1, 2, …, 5): 95% confidence region (thick dark line), mean of the stochastic response (thin dark line), and response of the nominal model (blue line).
Figures 5(a) and 6(a) show an increasing axial speed of the bit; it increases because of the constant axial force imposed at the bit (see Section 2.6). The behavior of the bit angular speed is similar for each Monte Carlo realization, and the values of this speed are between 0 and 350 rpm; Figs. 5(b) and 6(b) show that the upper confidence interval is close to the maximum response of the nominal model. Finally, Fig. 7 shows the random frequency spectrum, where a broad band spectrum is observed.
Let us consider now that μbit = 0.004 (the former bit-rock interaction used μbit = 0.04). This new value can be interpreted as the column facing a softer soil, easier to be drilled. In this case, the bit angular speed presents a different behavior as shown in Fig. 8, where three Monte Carlo samplings are drawn for the two different values of μbit. When μbit = 0.04 [former case, see Fig. 8(a)], the response has many oscillations and the stick-slip phenomenon is observed for all of the three Monte Carlo samplings. On the other hand, when μbit = 0.004 [see Fig. 8(b)], the response has less oscillations. In addition, for some combination of the parameters there is stick-slip, and for others there is not.
FIG. 8: Three Monte Carlo samplings of the bit angular speed for (a) μbit = 0.04 and (b) μbit = 0.004; δYi = 0.20 (i = 1, 2, …, 5).
If we define a stick-slip stability factor using the maximum and minimum values of the bit angular speed: 
=
(ωbmax − ωbmin) ⁄ (ωbmax + ωbmin), and a limit of
max = 0.80. Then, we can compute the probability of stickslip
to occur; i.e., the probability of to be greater than max. For example,
the probability of stick-slip to occur when μbit = 0.04 is 100%, and for μbit = 0.004, the probability of stick-slip is 50%.
If the probability of stick-slip occurrence was small, other techniques (such as subset simulations [24]) should be used to accelerate the stochastic
convergence, since many more samplings would be needed to achieve it.
5. CONCLUDING REMARKS
This paper has analyzed the nonlinear dynamics of a drill-string with uncertainties on the imposed speed and on the bit-rock parameters. The impact of these uncertainties on the dynamical response of the system has been investigated.
The continuum deterministic system is discretized by means of the finite-element method and a reduced-order model is constructed using the normal modes of the prestressed structure; only axial and torsional modes are included in the analysis. The results show that: (1) the initial prestressed configuration does not change significantly the natural frequencies and the normal modes of the structure; and (2) only two axial modes (one of which is a rigid body mode) and 10 torsional modes are necessary to represent the dynamics of the system, for the range of parameters analyzed.
Concerning the stochastic modeling, in the first step of the analysis, the parameters of the bit-rock interaction model, as well as the top speed, have been modeled with random variables and their probabilistic models have been constructed applying the maximum entropy principle. In the second step, the Monte Carlo method has been used to approximate the resulting stochastic differential equations.
The influence of each random variable on the system response is analyzed separately and also together, using different levels of uncertainty. Finally, the model has been used to compute the probability of occurrence of stickslip. Two different sets of parameters have been considered for the bit-rock interaction model. In the first set there has been a 100% probability of stick-slip, and in the second set there has been a 50% probability. From which, it can be concluded that, depending on the set of parameters of the system the sensitivity can be very high. In this sense, a careful identification procedure has to be done, such that the computational model can be used most effectively.
ACKNOWLEDGMENTS
The authors acknowledge the financial support of the Brazilian agencies CNPQ, CAPES, and FAPERJ.
REFERENCES
1. Spanos, P. D., Chevallier, A. M., Politis, N. P., and Payne, M. L., Oil and gas well drilling: A vibration perspective, Shock Vib. Dig., 35(2):85–103, 2003.
2. Spanos, P., Politis, N., Esteva, M., and Payne, M., Drillstring vibrations, in Advanced Drilling and Well Technology, pp. 117–156, Society of Petroleum Engineers, Richardson, TX, 2009.
3. Spanos, P. D., Chevallier, A. M., and Politis, N. P., Nonlinear stochastic drill-string vibrations, J. Vib. Acoust., 124(4):512–518, 2002.
4. Kotsonis, S. J. and Spanos, P. D., Chaotic and random whirling motion of drillstrings, J. Energy Resour. Technol., 119(4):217–222, 1997.
5. Ritto, T. G., Soize, C., and Sampaio, R., Nonlinear dynamics of a drill-string with uncertain model of the bit-rock interaction, Int. J. Non-Linear Mech., 44(8):865–876, 2009.
6. Ritto, T. G., Soize, C., and Sampaio, R., Probabilistic model identification of the bit-rock-interaction-model uncertainties in nonlinear dynamics of a drill-string, Mech. Res. Commun., 37(6):584–589, 2010.
7. Ritto, T. G., Soize, C., and Sampaio, R., Robust optimization of the rate of penetration of a drill-string using a stochastic nonlinear dynamical model, Comput. Mech., 45(5):415–427, 2010.
8. Ritto, T. G., Soize, C., and Sampaio, R., Stochastic dynamics of a drill-string with uncertain weight-on-hook, J. Braz. Soc. Mech. Sci. Eng., 32(3):250–258, 2010.
9. Soize, C., A nonparametric model of random uncertities for reduced matrix models in structural dynamics, Probab. Eng. Mech., 15:277–294, 2000.
10. Soize, C., Random matrix theory for modeling uncertainties in computational mechanics, Comput. Methods Appl. Mech. Eng., 194(12-16):1333–1366, 2005.
11. Yigit, A. and Christoforou, A., Coupled axial and transverse vibrations of oilwell drillstrings, J. Sound Vib., 195(4):617–627, 1996.
12. Yigit, A. and Christoforou, A., Coupled torsional and bending vibrations of drillstrings subject to impact with friction, J. Sound Vib., 215(1):167–181, 1998.
13. Christoforou, A. P. and Yigit, A. S., Fully coupled vibrations of actively controlled drillstrings, J. Sound Vib., 267:1029–1045, 2003.
14. Khulief, Y. A., Al-Sulaiman, F. A., and Bashmal, S., Vibration analysis of drillstrings with self excited stick-slip oscillations, J. Sound Vib., 299:540–558, 2007.
15. Trindade, M. A., Wolter, C., and Sampaio, R., Karhunen–Loève decomposition of coupled axial/bending of beams subjected to impacts, J. Sound Vib., 279:1015–1036, 2005.
16. Sampaio, R., Piovan, M. T., and Lozano, G. V., Coupled axial/torsional vibrations of drilling-strings by mean of nonlinear model, Mech. Res. Comun., 34(5-6):497–502, 2007.
17. Tucker, R. W. and Wang, C., An integrated model for drill-string dynamics, J. Sound Vib., 224(1):123–165, 1999.
18. Spanos, P. D., Sengupta, A. K., Cunningham, R. A., and Paslay, P. R., Modeling of roller cone bit lift-off dynamics in rotary drilling, J. Energy Res. Technol., 117(3):197–207, 1995.
19. Alamo, F. J. C., Dinâmica de estruturas unidimensionais esbeltas utilizando o contínuo de Cosserat, D.Sc. thesis, PUC-Rio, Rio de Janeiro, Brasil, 2006.
20. Jaynes, E., Probability Theory: The Logic of Science, vol. 1, Cambridge University Press, Cambridge, UK, 2003.
21. Shannon, C. E., A mathematical theory of communication, Bell Syst. Tech. J., 27:379–423 and 623–659, 1948.
22. Rubinstein, R. Y., Simulation and the Monte Carlo Method, Series in Probability and Statistics, 2nd ed., Wiley, Hoboken, NJ, 2007.
23. Idier, J., Bayesian Approach to Inverse Problems, Wiley-ISTE, Hoboken, NJ, 2008.
24. Au, S. K. and Beck, J. L., Estimation of small failure probabilities in high dimensions by subset simulation, J. Probab. Eng. Mech., 16(4):263–277, 2001.