# An Improved Sucker Rod Pumping System Model and Swabbing Parameters Optimized Design

Authors: Weicheng Li; Shimin Dong; Xiurong Sun
Journal: Mathematical Problems in Engineering
 (2018)

---

## Abstract

Considering the impact of fluid flowing into pump on sucker rod pumping system (SRPS) dynamic behaviors, an improved SRPS model with new boundary model is presented, which is a fluid-solid coupled model with the interactions among surface transmission, rod string longitudinal vibration, plunger motion, and fluid flow. A uniform algorithm is adopted instead of the mixed iteration algorithm for the surface transmission and downhole rod string vibration submodels, to reduce the difficulties of solving the entire SRPS model. The dynamic response comparison is executed between the improved model and the current model, and the results show that it will bring a calculation error on pump load and pump fullness if the progress of fluid flowing into the pump (PFFP) is ignored. Based on this improved model, a multitarget optimization model is proposed and the dynamic behavior of SRPS is improved with the optimized swabbing parameters.


---

## Body

## 1. Introduction

The SRPS is widely used in oil fields. It comprises three parts: surface transmission unit converting rotational motion into linear motion, sucker rod string as a joint between surface and downhole, and reciprocating pump exploiting the oil (see Figure1). The importance of predicting the dynamic responses of SRPS is to determine the operating situation and oil production [1]. For this equipment is usually set up in an open-air environment, and the test data device mounted on it is usually broken, especially the main working subsystem which is located nearly one kilometer or more downhole, making it difficult to be tested. Therefore, a more accurate SRPS simulation model should be established although the research in this field of study has been carried out widely. Due to the slender rod string moving upwards and downwards all time, an intense longitudinal vibration is produced. According to Figure 1, the SRPS model can be divided into rod string longitudinal vibration model, surface transmission model, and downhole pumping model. Commonly, this model is solved using the rod string longitudinal vibration equation as the foundation, surface and downhole model as boundary conditions. The most successful model of rod string longitudinal vibration is the Gibbs’s wave equation, and based on that, the models are studied specifically on the enhancement of the surface and downhole boundaries with different operating conditions [2–4]. The surface boundary condition that includes the motor speed variations has been extensively used and is more applicable in practice [5, 6] whereas the downhole boundary condition has been continuously improved. However, further study is still required due to the inconsistent alteration and complication encountered in downhole operation.Figure 1
Sucker rod pumping system.The downhole boundary condition actually is a model describing the pump operation. The pump operating mode can be divided into upstroke and downstroke. During upstroke, as the plunger moves upwards, the pump pressure will decrease and the fluid will not be sucked into the pump until this pressure drop down to pump inlet pressure. For the downstroke, the pressure will increase with the plunger goes down and the fluid will begin to be drained out when the pressure equals the pump outlet pressure. Based on the SRPS’s operating state, the first universal downhole boundary model was divided into four phases with a vague formulation [2]. In an ideal condition, during the fully loaded upstroke movement, the pump load was set equivalent to the fluid load; pump load was set to zero for the unloaded downstroke movement [7, 8]. As the pump operation phases are highly related to pump pressure, it was revised with an explicit formation deduced by the interaction between pump outlet pressure and the pressure in the pump barrel [9–14], as follows: (1)FPLt=Appd-p-ArpdwhereFPL is the pump load, N;Ap andAr are the cross section area of plunger and rod string respectively, m2;pd is the pump outlet pressure, pa;p is the pump pressure, pa.For thepd is always considered as a constant pressure, the key of this research is to establish an accurate pump pressure model which consists of four phases as shown in Figure 2. With considering of the gas, this downhole boundary model is improved as shown in (2). When the pump is at phases 1 and 3, the pressure variation obeys the rules of gas state equation. As for in phases 2 and 4, the pump pressure is taken as theps andpd separately.(2)p=LogLog+up-∫0tsq/ApQpdphase1p=psphase2p=LgLg-Ls+up-∫tuttq/ApQpsphase3p=pdphase4where Log and Lg are the gas column length when plunger is arriving at bottom dead center and top dead center, respectively, m;up is the plunger displacement, m;Ls is the pump stroke displacement, m;Lp is the plunger length, m;μ is the fluid dynamic viscosity, pa·s;δ is the clearance between plunger and pump barrel, m;Dd is the pump diameter, m;ts andtt are the open time of standing valve and travelling valve, respectively, s;tu is the upstroke time;q is the liquid instantaneous leakage volume.Figure 2
Four phases of pump operation.In this formula, the principle that the gas/oil ratio of clearance volume (the space volume when plunger arrives at bottom dead center) equals the gas/oil ratio of pump inlet is applied. However, at this time, the pump pressure should be the same as the pump outlet pressure. Thus, it is revised with the pump outlet pressure and this improved model has been widely used until now. [15–18]. However, the current downhole boundary model still exists some shortcomings need to be improved. Due to ignoring the progress of PFFP, it is established with the hypothesis of regarding the pump is filled with the fluid, whose gas/liquid ratio always equals that at pump intake, as well as the one in the clearance volume. This assumption is not applicable for the oil well with insufficient oil well deliverability (OWD) which will result in the pumping fluid being unable to keep pace with the plunger and cause an incomplete fullness. Meanwhile, a pump load calculation error will be produced with keeping the pump pressure as constant when fluid flows into the pump.The SRPS model lays the foundation for optimizing the swabbing parameters to improve the system operation status, except for predicting and evaluating its dynamic response. Miska et al., 1997 [19], propose a computer-aided optimization method relying on a simple linear algebraic system model, so as to minimize the energy consumption. Firu et al., 2003 [20], present an improve optimization criterion including eight operational parameters to achieve the maximum system efficiency. Liu and Qi, 2011 [21], apply fluid flow characteristics in coalbed methane reservoirs to estimate the production capacity and build the system efficiency optimization model combined with SRPS performance. The above optimization models are built with a simplified pump load; then an improved system efficiency optimization model is built jointing with formula (2) [22, 23]. However, the current optimization models ignore the effect of swabbing parameters on pump fullness, in that the SRPS model’s restriction. Besides that, single target optimization cannot make an accurate and comprehensive presentation for SRPS whose pumping progress is complex, multicomponent and interactive.In this paper, firstly, an improved SRPS model is presented with the new downhole boundary model, considering the movement of fluid flowing into pump with gas instantaneous dissolution and evolution. Secondly, for the new downhole boundary model, a nonlinear fluid-solid coupled model will increase the complexity of this improved SRPS model. A unified numerical algorithm is applied on the whole model to decrease the calculation time, instead of traditional mixed iteration algorithm. Thirdly, a multiobjective optimization model is proposed based on the improved SRPS model. Fourthly, the surface dynamometer card is collected to verify the improved model’s accuracy. Fifthly, the dynamic response comparison on the current SRPS model and improved SRPS model is executed. Finally, the optimization program is applied on a test well, and the results are given before and after optimization.


## 2. Integrated Simulation Model

The improved SRPS model is a nonlinear fluid-solid coupled model consisted of drive, transmission and load. The motor is used to generate power and its rotation motion is transformed into linear motion of suspension point. This linear motion is passed on to the plunger through rod string longitudinal vibration; meanwhile, plunger motion has an influence on pumping fluid by regulating pump pressure. On the contrary, fluid motion affects pump pressure which is also the main factor in determining pump load. This pump load and other system loads are transferred to motor output shaft by surface and downhole transmission, causing impact on the motor motion. The above coupled process is depicted in Figure3. In this paper, the SRPS model is subdivided into rod string longitudinal vibration model, corresponding surface boundary and downhole boundary model, so as to have a clear description.Figure 3
SRPS coupled model sketch.


### 2.1. Sucker Rod Longitudinal Vibration Model

While the SRPS is operating, the motion of arbitrary micro element of the slender rod string produced intense longitudinal vibration can be composed of the following two parts: (1) following the movement of suspension point and (2) relative movement to the suspension point. Referring to Figure4, the coordinate system is built based on the top dead center of horse head as the origin. Then the stress balance terms of rod string element and top and bottom boundary conditions are used to establish the sucker rod longitudinal vibration model, with the assumption of ignoring the deformation of tube and the vibration of liquid column in a vertical well.(3)∂2u∂t2-Erρr∂2u∂x2+μρrAr∂u∂t=-d2uatdt2-μρrArduatdt+gux,tx=0=uatErAr∂ux,t∂xx=Lr=FPLtwhereu is the displacement of the rod string at arbitrary depth and time, m;Er andρr are the elasticity modulus and density of rod string, respectively, pa and kg/m3;ua is the suspension point displacement, m.Figure 4
Mechanical model of rod string longitudinal vibration.


### 2.2. Surface Boundary Model

The surface boundary condition refers to the motion of suspension point as in (3). With the assumption that the surface transmission mechanism consists of belt-reducing gear and four-bar linkage as a single degree of freedom system, it can be represented as a function of motor angle. Later, the surface apparatus motion model is set up using Lagrange equation, taking motor output shaft as equivalent component, its motion as generalized coordinates and set 12 o’clock as a reference direction (Figure 5).(4)Ieθ¨m+12θ˙m2dIedθm=Med-Mefθmt=0=0θ˙mθ=0=ω0Figure 5
Surface apparatus equivalent motion model.With the equivalent rotating inertia(5)Ie=∑Ijωjθ˙m2+∑mjvjθ˙m2The semiempirical formulation of motor driving torqueMed deduced from the speed/torque characteristic is adopted from Wu et al. [24].(6)Med=2λkMHεcωnωn-θ˙mεc2ωn2+ωn-θ˙m2whereλk is the motor overload coefficient;εc is the motor critical-slip ratio;ωn is the motor synchronous angular velocity, rad/s.(7)MH=9550PHnHεc=εrλk+λk2-1ωn=2πfzwherez is the motor pole pairs;f is the power frequency, Hz;PH is the motor rated power, kW;nH is the motor rated speed, min−1;εr is the motor rated-slip ratio.The equivalent resistance torque is derived in accordance with the force balance of four-bar linkage (Figure6), written as(8)Mef=1ibgεbSTfFRL-Fbεpk1-Mcsin⁡θ-θτFigure 6
Four-bar linkage.with the angle of crank(9)θ=θmibgwhereibg is the transmission ratio of the belt-reducing gear system;εp is the transmission efficiency from crank to suspension point;STf is torque factor, m;FRL is the suspension point load, N;Fb is the structural unbalance weight, N;Mc is the maximum balancing torque of the crank N·m.The displacement of suspension point can be expressed as motor angle, based on the geometry relation shown in Figure6 [6]. The surface boundary motion model is(10)uat=arccos⁡LC2+LK2-LR+LP22LCLK-arccos⁡LC2+LL2-LP22LCLL-arcsin⁡LRLLsin⁡2π-θ0-θmibg+arcsin⁡LILK


### 2.3. Downhole Boundary Model

This new model is an improved method from the current boundary models considering the interaction among pump pressure, fluid flow into pump, gas dissolution and evolution, pressure drop due to fluid gravitational potential energy, inertia head loss, friction head loss, and local head loss. The new boundary model consists of pressure variation model and fluid flow into pump model. Some assumptions are made for building this model:(1)
Bottom dead point of plunger as the origin position.(2)
The pump inlet pressureps and outlet pressurepd are kept constant for one cycle of SRPS operation.Pressure Variation Model. The gas state equations are built as follows ( (11) and (12)) when the displacement of plunger isup orup+dup:(11)pup+Log-Lf-LkAp=ZNQT(12)p+dpup+dup+Log-Lf-Lk-dLf-dLkAp=Z+∂Z∂pdpN+dNQTwhereN is the gas molar number, mol;Q is the natural gas constant, J/(mol·K);T is the temperature in pump barrel, °C; Z is the natural gas compressibility factor;Lk is the equivalent length of pump leakage, m;Lf is the liquid level in the pump, m.Dividing (11) by (12) and ignoring the second order small quantities, then the variation of pump pressure with plunger displacement is obtained:(13)dpdup=1/NdN/dup+1/up+Log-Lf-LkdLf+Lk/dup-1p-1-1/ZdZ/dpThe above equation is converted into a time varying function with the purpose of facilitating it with the rod string longitudinal vibration equation.(14)dpdt=1/NdN/dt1/vp+1/up+Log-Lf-LkdLf+Lk/dt1/vp-1p-1-1/ZdZ/dpvpwherevp is the plunger velocity, m/s.From (14),dNis composed of the following sections.(1) The molar number of gas is released from or dissolved into oil as pressure varies(15)dNr=-ZTpstZstTstpα·dp103ApLf+Lo-Log1-εwBoVmolwherepst is the standard pressure, pa;Zst is the standard natural gas compressibility factor;Tst is the standard temperature, °C; α is the solubility coefficient of natural gas m3/ (m3·pa);εw is the water content;Vmol is the molar volume, m3;Bo is the crude oil volume factor.(2) The molar number of gas with fluid being sucked into the pump corresponding to the pump pressure is(16)dNs=103ApRJdLfVmolwhereRJ is the transient gas liquid ratio in the pump, m3/m3.(3) The molar number of gas with leaked fluid at the pump pressure is(17)dNv=103ApRJdLvVmolwhere(18)RJ=1-εwRp-RspstZTTstZstpRs=αp-pstwhereRp is the production gas oil ratio, m3/m3;Rs is the gas oil ratio at pump intake, m3/m3.Fluid Flowing into Pump Model. Figure 7 shows the pump operation when the gas-liquid flow is drawn into the pump. Sectiona-a indicates the cross section of standing valve hole. Supposing the sectionf-fis the liquid level at arbitrary time, neglect the process of gas bubbling from fluid by considering only the pump gas to occupy the upper space off-f. The one-dimensional unsteady flow equation based on Bernoulli equation is established to describe the flow state of pumping liquid. In the Bernoulli equation, the inertia head loss, friction head loss, and local head loss are all considered. (19)psρg+va22g=pρg+vf22g+Lo-Log+Lf+Lv+ha+hv+hfwhere(20)ha=1gLo-Log+Lf+LvAfApdvfdthf=64μvfLo-Log+Lf2Dd2ρghv=1gε2ApAf2va22va=ApAfvfwhereε is the flow coefficient of standing valve.Figure 7
Process of fluid flows into the pump.Then (19) is converted into differential forms.(21)L˙f=vfv˙f=L0-VogAp+Lf+LvAfAp-1psρ-pρ+vf22Ap2Af2-Ap2ε2Af2-64μL0-Vog/Ap+LfDd2vfρ-1-gL0-VogAp+Lf+LvNormally, the fluid in the barrel, no matter in which form, will all be discharged except the one in dead space. Therefore, the fluid flowing out of the model with opening traveling valve is not proposed. The pump outlet pressure generated by several kilometers liquid column is large compared to the hydraulic loss, causing the ignorance of hydraulic loss when fluid flows out. Based on the statement mentioned above, the new downhole boundary model is expressed as follows and also is divided into four phases corresponding to the phases shown in Figure2:(22)Status1:beforeopeningthestandingvalvedpdt=q/Ap·vp1/up-Lv+Log+103Ap·RJ/N·Vmol-1/up-Lv+Logp-1-1/Z·dZ/dp+103α/N·ZTpst/Tstp·Lv+Lo-Log/spVmol·Ap·1-εw·vpStatus2:afteropeningthestandingvalvedpdt=vf·Ap+q/Ap·vp1/up-Lf-Lv+Log+103Ap·RJ/N·Vmol-1/up-Lf-Lv+Logp-1-1/Z·dZ/dp+103α/N·ZTpst/Tstp·Lv+Lf+Lo-Log/spVmol·Ap·1-εw·vpStatus3:beforeopeningthetravelingvalvedpdt=q/Ap·vp1/up-Lf-Lv+Log+103Ap·RJ/N·Vmol-1/up-Lf-Lv+Logp-1-1/Z·dZ/dp+103α/N·ZTpst/Tstp·Lv+Lf+Lo-Log/spVmol·Ap·1-εw·vpStatus4:afteropenginthetravelingvalvedpdt=0where(23)Vmol=22.4ZTpstTstp


## 2.1. Sucker Rod Longitudinal Vibration Model

While the SRPS is operating, the motion of arbitrary micro element of the slender rod string produced intense longitudinal vibration can be composed of the following two parts: (1) following the movement of suspension point and (2) relative movement to the suspension point. Referring to Figure4, the coordinate system is built based on the top dead center of horse head as the origin. Then the stress balance terms of rod string element and top and bottom boundary conditions are used to establish the sucker rod longitudinal vibration model, with the assumption of ignoring the deformation of tube and the vibration of liquid column in a vertical well.(3)∂2u∂t2-Erρr∂2u∂x2+μρrAr∂u∂t=-d2uatdt2-μρrArduatdt+gux,tx=0=uatErAr∂ux,t∂xx=Lr=FPLtwhereu is the displacement of the rod string at arbitrary depth and time, m;Er andρr are the elasticity modulus and density of rod string, respectively, pa and kg/m3;ua is the suspension point displacement, m.Figure 4
Mechanical model of rod string longitudinal vibration.


## 2.2. Surface Boundary Model

The surface boundary condition refers to the motion of suspension point as in (3). With the assumption that the surface transmission mechanism consists of belt-reducing gear and four-bar linkage as a single degree of freedom system, it can be represented as a function of motor angle. Later, the surface apparatus motion model is set up using Lagrange equation, taking motor output shaft as equivalent component, its motion as generalized coordinates and set 12 o’clock as a reference direction (Figure 5).(4)Ieθ¨m+12θ˙m2dIedθm=Med-Mefθmt=0=0θ˙mθ=0=ω0Figure 5
Surface apparatus equivalent motion model.With the equivalent rotating inertia(5)Ie=∑Ijωjθ˙m2+∑mjvjθ˙m2The semiempirical formulation of motor driving torqueMed deduced from the speed/torque characteristic is adopted from Wu et al. [24].(6)Med=2λkMHεcωnωn-θ˙mεc2ωn2+ωn-θ˙m2whereλk is the motor overload coefficient;εc is the motor critical-slip ratio;ωn is the motor synchronous angular velocity, rad/s.(7)MH=9550PHnHεc=εrλk+λk2-1ωn=2πfzwherez is the motor pole pairs;f is the power frequency, Hz;PH is the motor rated power, kW;nH is the motor rated speed, min−1;εr is the motor rated-slip ratio.The equivalent resistance torque is derived in accordance with the force balance of four-bar linkage (Figure6), written as(8)Mef=1ibgεbSTfFRL-Fbεpk1-Mcsin⁡θ-θτFigure 6
Four-bar linkage.with the angle of crank(9)θ=θmibgwhereibg is the transmission ratio of the belt-reducing gear system;εp is the transmission efficiency from crank to suspension point;STf is torque factor, m;FRL is the suspension point load, N;Fb is the structural unbalance weight, N;Mc is the maximum balancing torque of the crank N·m.The displacement of suspension point can be expressed as motor angle, based on the geometry relation shown in Figure6 [6]. The surface boundary motion model is(10)uat=arccos⁡LC2+LK2-LR+LP22LCLK-arccos⁡LC2+LL2-LP22LCLL-arcsin⁡LRLLsin⁡2π-θ0-θmibg+arcsin⁡LILK


## 2.3. Downhole Boundary Model

This new model is an improved method from the current boundary models considering the interaction among pump pressure, fluid flow into pump, gas dissolution and evolution, pressure drop due to fluid gravitational potential energy, inertia head loss, friction head loss, and local head loss. The new boundary model consists of pressure variation model and fluid flow into pump model. Some assumptions are made for building this model:(1)
Bottom dead point of plunger as the origin position.(2)
The pump inlet pressureps and outlet pressurepd are kept constant for one cycle of SRPS operation.Pressure Variation Model. The gas state equations are built as follows ( (11) and (12)) when the displacement of plunger isup orup+dup:(11)pup+Log-Lf-LkAp=ZNQT(12)p+dpup+dup+Log-Lf-Lk-dLf-dLkAp=Z+∂Z∂pdpN+dNQTwhereN is the gas molar number, mol;Q is the natural gas constant, J/(mol·K);T is the temperature in pump barrel, °C; Z is the natural gas compressibility factor;Lk is the equivalent length of pump leakage, m;Lf is the liquid level in the pump, m.Dividing (11) by (12) and ignoring the second order small quantities, then the variation of pump pressure with plunger displacement is obtained:(13)dpdup=1/NdN/dup+1/up+Log-Lf-LkdLf+Lk/dup-1p-1-1/ZdZ/dpThe above equation is converted into a time varying function with the purpose of facilitating it with the rod string longitudinal vibration equation.(14)dpdt=1/NdN/dt1/vp+1/up+Log-Lf-LkdLf+Lk/dt1/vp-1p-1-1/ZdZ/dpvpwherevp is the plunger velocity, m/s.From (14),dNis composed of the following sections.(1) The molar number of gas is released from or dissolved into oil as pressure varies(15)dNr=-ZTpstZstTstpα·dp103ApLf+Lo-Log1-εwBoVmolwherepst is the standard pressure, pa;Zst is the standard natural gas compressibility factor;Tst is the standard temperature, °C; α is the solubility coefficient of natural gas m3/ (m3·pa);εw is the water content;Vmol is the molar volume, m3;Bo is the crude oil volume factor.(2) The molar number of gas with fluid being sucked into the pump corresponding to the pump pressure is(16)dNs=103ApRJdLfVmolwhereRJ is the transient gas liquid ratio in the pump, m3/m3.(3) The molar number of gas with leaked fluid at the pump pressure is(17)dNv=103ApRJdLvVmolwhere(18)RJ=1-εwRp-RspstZTTstZstpRs=αp-pstwhereRp is the production gas oil ratio, m3/m3;Rs is the gas oil ratio at pump intake, m3/m3.Fluid Flowing into Pump Model. Figure 7 shows the pump operation when the gas-liquid flow is drawn into the pump. Sectiona-a indicates the cross section of standing valve hole. Supposing the sectionf-fis the liquid level at arbitrary time, neglect the process of gas bubbling from fluid by considering only the pump gas to occupy the upper space off-f. The one-dimensional unsteady flow equation based on Bernoulli equation is established to describe the flow state of pumping liquid. In the Bernoulli equation, the inertia head loss, friction head loss, and local head loss are all considered. (19)psρg+va22g=pρg+vf22g+Lo-Log+Lf+Lv+ha+hv+hfwhere(20)ha=1gLo-Log+Lf+LvAfApdvfdthf=64μvfLo-Log+Lf2Dd2ρghv=1gε2ApAf2va22va=ApAfvfwhereε is the flow coefficient of standing valve.Figure 7
Process of fluid flows into the pump.Then (19) is converted into differential forms.(21)L˙f=vfv˙f=L0-VogAp+Lf+LvAfAp-1psρ-pρ+vf22Ap2Af2-Ap2ε2Af2-64μL0-Vog/Ap+LfDd2vfρ-1-gL0-VogAp+Lf+LvNormally, the fluid in the barrel, no matter in which form, will all be discharged except the one in dead space. Therefore, the fluid flowing out of the model with opening traveling valve is not proposed. The pump outlet pressure generated by several kilometers liquid column is large compared to the hydraulic loss, causing the ignorance of hydraulic loss when fluid flows out. Based on the statement mentioned above, the new downhole boundary model is expressed as follows and also is divided into four phases corresponding to the phases shown in Figure2:(22)Status1:beforeopeningthestandingvalvedpdt=q/Ap·vp1/up-Lv+Log+103Ap·RJ/N·Vmol-1/up-Lv+Logp-1-1/Z·dZ/dp+103α/N·ZTpst/Tstp·Lv+Lo-Log/spVmol·Ap·1-εw·vpStatus2:afteropeningthestandingvalvedpdt=vf·Ap+q/Ap·vp1/up-Lf-Lv+Log+103Ap·RJ/N·Vmol-1/up-Lf-Lv+Logp-1-1/Z·dZ/dp+103α/N·ZTpst/Tstp·Lv+Lf+Lo-Log/spVmol·Ap·1-εw·vpStatus3:beforeopeningthetravelingvalvedpdt=q/Ap·vp1/up-Lf-Lv+Log+103Ap·RJ/N·Vmol-1/up-Lf-Lv+Logp-1-1/Z·dZ/dp+103α/N·ZTpst/Tstp·Lv+Lf+Lo-Log/spVmol·Ap·1-εw·vpStatus4:afteropenginthetravelingvalvedpdt=0where(23)Vmol=22.4ZTpstTstp


## 3. Integrated Numerical Algorithm

For the current SRPS model, the surface transmission model is solved by numerical integration method [25]. While the finite difference method is used to describe the rod string longitudinal vibration [2, 26, 27]. Since there is no fixed solution between surface and downhole model, it needs cyclic iteration to handle the whole model, which is time consuming [25]. Considering the interaction within plunger motion, fluid flow into pump and pump pressure, the improved SRPS model has a higher nonlinear degree and thus increases the difficulty of solving it. In order to reduce the solving time, the wave equation of rod string longitudinal vibration is converted into ordinary differential equations by modal superposition method. Then the surface and downhole boundary model are now in the form of ordinary differential equations at which the whole model can be solved by Runge-Kutta method directly. Besides that, solving the rod string longitudinal vibration equation is also solving the force vibration response of rod string, where (3) is converted into the following form:(24)ρrAr∂2u∂t2-ρrEr∂2u∂x2+μ∂u∂t=ρrArfx,twhere(25)fx,t=-d2uatdt2-μρrArduatdt+gAnd its solution can be expressed in the function of time and space by separating the variables:(26)ux,t=∑i=1∞ΦixqitTherefore, (24) can be expressed as(27)∑i=1∞ArρrΦixq¨it-∑i=1∞ErArd2Φixdx2qit+∑i=1∞μΦixq˙it=Arρrfx,tMultiply the above equation byΦi(x), and it is integrated along the rod length. Then the following is derived:(28)Arρrq¨it∫0LrΦi2xdx-ErArqit∫0LrΦixΦi′′xdx+μq˙it∫0LrΦi2xdx=Arρr∫0LrΦixfx,tdxIts mode shapes and natural frequencies are(29)Φix=2ρrArsin⁡2i-1π2Lrxpni=2i-1πEr/ρr2LrThen (28) can be simplified as follows using mode shape orthogonality:(30)q¨it+μArρrq˙it+2i-1π2Lr2Erρrqit=Fiwhere(31)Fi=d2uatdt2+μduatdt22ArρrLr2i-1π+FPLt2ArρrLrsin⁡2i-1π2Letith-order forced vibration displacement response and velocity response bexi1 andxi2, respectively. Then (30) can be expressed as the following form:(32)x˙i1t=xi2tx˙i2t=-μArρrxi2t-2i-1π2Lr2Erρrxi1t+FiThen suspension point load is derived.(33)FRLt=ErAr∂u0,t∂x+Fr=ErAr2ρrArLr∑i=1N2i-1π2Lrxi1t+FrThe displacement and velocity of plunger are(34)upt=uat-uL,t=uat-2ArρrLr∑i=1i=2∗j+1xi1t-∑i=1i=2∗jxi1tvpt=u˙at-∂uL,t∂t=u˙at-2ArρrLr∑i=1i=2∗j+1xi2t-∑i=1i=2∗jxi2tThen the integrated numerical model is given as follows:(35)
θ
˙
m
=
ω
m
ω
˙
m
=
1
I
e
2
λ
k
M
H
ε
c
ω
n
ω
n
-
θ
˙
m
ε
c
2
ω
n
2
+
ω
n
-
θ
˙
m
2
-
1
i
b
g
η
b
g
S
T
f
F
R
L
-
F
b
η
C
L
k
1
-
M
c
sin
⁡
θ
-
θ
τ
-
1
2
ω
m
2
d
I
e
d
θ
m
u
a
t
=
arccos
⁡
L
C
2
+
L
K
2
-
L
R
+
L
P
2
2
L
C
L
K
-
arccos
⁡
L
C
2
+
L
L
2
-
L
P
2
2
L
C
L
L
-
arcsin
⁡
L
R
L
L
sin
⁡
2
π
-
θ
0
-
θ
m
i
b
g
+
arcsin
⁡
L
I
L
K
x
11
t
=
x
12
t
x
˙
12
t
=
-
μ
A
r
ρ
r
x
12
t
-
π
2
L
r
2
E
r
ρ
r
x
11
t
+
d
2
u
a
t
d
t
2
+
μ
d
u
a
t
d
t
2
2
A
r
ρ
r
L
r
π
+
A
p
p
d
-
p
-
A
r
p
d
+
F
f
2
A
r
ρ
r
L
r
x
˙
21
t
=
x
22
t
x
˙
22
t
=
-
μ
A
r
ρ
r
x
22
t
-
3
π
2
L
r
2
E
r
ρ
r
x
21
t
+
d
2
u
a
t
d
t
2
+
μ
d
u
a
t
d
t
2
2
A
r
ρ
r
L
r
3
π
-
A
p
p
d
-
p
-
A
r
p
d
+
F
f
2
A
r
ρ
r
L
r
⋮
x
˙
i
1
t
=
x
i
2
t
x
˙
i
2
t
=
-
μ
A
r
ρ
r
x
i
2
t
-
2
i
-
1
π
2
L
r
2
E
r
ρ
r
x
i
1
t
+
d
2
u
a
t
d
t
2
+
μ
d
u
a
t
d
t
2
2
A
r
ρ
r
L
r
2
i
-
1
π
+
A
p
p
d
-
p
-
A
r
p
d
+
F
f
2
A
r
ρ
r
L
r
sin
⁡
2
i
-
1
π
2
F
R
L
t
=
E
r
A
r
2
ρ
r
A
r
L
r
∑
i
=
1
N
2
i
-
1
π
2
L
r
x
i
1
t
+
F
r
u
p
t
=
u
a
t
-
2
A
r
ρ
r
L
r
∑
i
=
1
i
=
2
∗
j
+
1
x
i
1
t
-
∑
i
=
1
i
=
2
∗
j
x
i
1
t
v
p
t
=
u
˙
a
t
-
2
A
r
ρ
r
L
r
∑
i
=
1
i
=
2
∗
j
+
1
x
i
2
t
-
∑
i
=
1
i
=
2
∗
j
x
i
2
t
S
t
a
t
u
s
1:
b
e
f
o
r
e
o
p
e
n
i
n
g
t
h
e
s
t
a
n
d
i
n
g
v
a
l
v
e
d
p
d
t
=
q
/
A
p
·
v
p
1
/
u
p
-
L
v
+
L
o
g
+
1
0
3
A
p
·
R
J
/
N
·
V
m
o
l
-
1
/
u
p
-
L
v
+
L
o
g
p
-
1
-
1
/
Z
·
d
Z
/
d
p
+
1
0
3
α
/
N
·
Z
T
p
s
t
/
T
s
t
p
·
L
v
+
L
o
-
L
o
g
/
s
p
V
m
o
l
·
A
p
·
1
-
ε
w
·
v
p
S
t
a
t
u
s
2:
a
f
t
e
r
o
p
e
n
i
n
g
t
h
e
s
t
a
n
d
i
n
g
v
a
l
v
e
d
p
d
t
=
v
f
·
A
p
+
q
/
A
p
·
v
p
1
/
u
p
-
L
f
-
L
v
+
L
o
g
+
1
0
3
A
p
·
R
J
/
N
·
V
m
o
l
-
1
/
u
p
-
L
f
-
L
v
+
L
o
g
p
-
1
-
1
/
Z
·
d
Z
/
d
p
+
1
0
3
α
/
N
·
Z
T
p
s
t
/
T
s
t
p
·
L
v
+
L
f
+
L
o
-
L
o
g
/
s
p
V
m
o
l
·
A
p
·
1
-
ε
w
·
v
p
S
t
a
t
u
s
3:
b
e
f
o
r
e
o
p
e
n
i
n
g
t
h
e
t
r
a
v
e
l
i
n
g
v
a
l
v
e
d
p
d
t
=
q
/
A
p
·
v
p
1
/
u
p
-
L
f
-
L
v
+
L
o
g
+
1
0
3
A
p
·
R
J
/
N
·
V
m
o
l
-
1
/
u
p
-
L
f
-
L
v
+
L
o
g
p
-
1
-
1
/
Z
·
d
Z
/
d
p
+
1
0
3
α
/
N
·
Z
T
p
s
t
/
T
s
t
p
·
L
v
+
L
f
+
L
o
-
L
o
g
/
s
p
V
m
o
l
·
A
p
·
1
-
ε
w
·
v
p
S
t
a
t
u
s
4:
a
f
t
e
r
o
p
e
n
i
n
g
t
h
e
t
r
a
v
e
l
i
n
g
v
a
l
v
e
d
p
d
t
=
0
d
L
f
d
t
=
v
f
d
v
f
d
t
=
L
0
-
V
o
g
A
p
+
L
f
+
L
v
A
f
A
p
-
1
p
s
ρ
-
p
ρ
+
v
f
2
2
A
p
2
A
f
2
-
A
p
2
ε
2
A
f
2
-
64
μ
L
0
-
V
o
g
/
A
p
+
L
f
D
2
v
f
ρ
-
1
-
g
L
0
-
V
o
g
A
p
+
L
f
+
L
v


## 4. Optimization Model

### 4.1. Optimization Goal

As mostly developed oil-field moves into the mid and late stage and the OWD begins to decline gradually, energy-saving, production-increasing, and reducing load variation as much as possible are particularly important. Then we take pump fullnessepf, suspension point load amplitudeFRLA, crank torque standard deviationMcsd, and motor input power average P-m as the optimization goal to build a multitarget model. Suspension point load and crank torque can be solved directly by (35). The motor input power and pump fullness calculation formula are deduced as follows:(36)Pm=Medθ˙m+P0+1ηH-1PH-P0Medθ˙mPH2(37)epf=Lf-∫0tuqdt-Lo-LogupuwhereP0 is the no-load power of motor, kW;PH is the motor rated power, kW;ηH is the motor rated efficiency, kW;upu is the plunger displacement when it is arriving at dead top center.So the objective function is(38)ΩepfX,FRLAX,McsdX,P-mX=K1epfX+K2FRLAX+K3McsdX+K4P-mXwhereK1,K2,K3, andK4 are the weight coefficients.


### 4.2. Design Variables

The objective function can be expressed as the function of the swabbing parameters when the SRPS type and oil well basic parameters are confirmed. In this paper, the swabbing parameters denote strokeS, stroke frequencyns, pump diameterDd, pump depthLpd, crank balance radiusrc, and rod string combination (thejth rod string diameter and length aredj andLj,, respectively, k=1,2,…,m). Then the design variables are shown as follows:(39)X=S,ns,Dd,Lpd,dj,Lj;j=1,2,⋯m


### 4.3. Constraints

(a) Crank Balance Degree. Crank balance degree indicates the load fluctuation to a certain degree and it needs to be kept at a high value: (40)0.95≤MckuMckd≤1whereMcku andMckd are the maximum crank torque when plunger is at upstroke and downstroke, respectively, N·m.(b) Ground Device Carrying Capacity. The suspension point load, crank torque and motor torque, at any time [0,T] do not go beyond the allowable range {max⁡FRL,min⁡FRL;max⁡Mc,min⁡Mc;max⁡Med, min⁡(Med)} for the given type: (41)min⁡FRL≤FRL≤max⁡FRLmin⁡Mck≤Mck≤max⁡Mckmin⁡Med≤Med≤max⁡Med(c) Rod String Strength. The maximum and minimum stress of any pointx along the rod string does not exceed permissible stress range in one cycle: (42)σmin≤σrs≤σmax(d) Swabbing Parameters. Each oil well has different limit on the swabbing parameters in accordance with the device type and actual operation; hence their allowable variation ranges are adjustable.


### 4.4. Optimization Algorithm

In summary, this multivariable optimization model is established with nonlinear restriction and nonlinear objective function. So as to seek the best results, the genetic algorithm is applied to solve it.


## 4.1. Optimization Goal

As mostly developed oil-field moves into the mid and late stage and the OWD begins to decline gradually, energy-saving, production-increasing, and reducing load variation as much as possible are particularly important. Then we take pump fullnessepf, suspension point load amplitudeFRLA, crank torque standard deviationMcsd, and motor input power average P-m as the optimization goal to build a multitarget model. Suspension point load and crank torque can be solved directly by (35). The motor input power and pump fullness calculation formula are deduced as follows:(36)Pm=Medθ˙m+P0+1ηH-1PH-P0Medθ˙mPH2(37)epf=Lf-∫0tuqdt-Lo-LogupuwhereP0 is the no-load power of motor, kW;PH is the motor rated power, kW;ηH is the motor rated efficiency, kW;upu is the plunger displacement when it is arriving at dead top center.So the objective function is(38)ΩepfX,FRLAX,McsdX,P-mX=K1epfX+K2FRLAX+K3McsdX+K4P-mXwhereK1,K2,K3, andK4 are the weight coefficients.


## 4.2. Design Variables

The objective function can be expressed as the function of the swabbing parameters when the SRPS type and oil well basic parameters are confirmed. In this paper, the swabbing parameters denote strokeS, stroke frequencyns, pump diameterDd, pump depthLpd, crank balance radiusrc, and rod string combination (thejth rod string diameter and length aredj andLj,, respectively, k=1,2,…,m). Then the design variables are shown as follows:(39)X=S,ns,Dd,Lpd,dj,Lj;j=1,2,⋯m


## 4.3. Constraints

(a) Crank Balance Degree. Crank balance degree indicates the load fluctuation to a certain degree and it needs to be kept at a high value: (40)0.95≤MckuMckd≤1whereMcku andMckd are the maximum crank torque when plunger is at upstroke and downstroke, respectively, N·m.(b) Ground Device Carrying Capacity. The suspension point load, crank torque and motor torque, at any time [0,T] do not go beyond the allowable range {max⁡FRL,min⁡FRL;max⁡Mc,min⁡Mc;max⁡Med, min⁡(Med)} for the given type: (41)min⁡FRL≤FRL≤max⁡FRLmin⁡Mck≤Mck≤max⁡Mckmin⁡Med≤Med≤max⁡Med(c) Rod String Strength. The maximum and minimum stress of any pointx along the rod string does not exceed permissible stress range in one cycle: (42)σmin≤σrs≤σmax(d) Swabbing Parameters. Each oil well has different limit on the swabbing parameters in accordance with the device type and actual operation; hence their allowable variation ranges are adjustable.


## 4.4. Optimization Algorithm

In summary, this multivariable optimization model is established with nonlinear restriction and nonlinear objective function. So as to seek the best results, the genetic algorithm is applied to solve it.


## 5. Test and Verification

Surface dynamometer card is a closed graph recording polished rod loads versus rod displacement over a SRPS cycle, which is generally collected and taken as an index to estimate the operation of SRPS. In this paper, based on the dynamometer sensor shown in Figure8, four test wells are used to validate the improved SRPS model. The oil well parameters are listed in Table 1, and the simulation and field test results are given in Figure 9. According to the plotted curves, the simulation results are basically consistent with that in measured, and the improved SRPS model is accurate enough to be proposed for engineering practices.Table 1
Oil well basic parameters.
Well #1 #2 #3 #4 Motor type YD280S-8 Y250M-6 Y250M-6 YD280S-6 Pumping unit type CYJ10-3-53HB CYJ10-3-53HB CYJ10-3-53HB CYJ10-3-53HB Stroke length (m) 3 3 3 3 Stroke frequency (min−1) 3.5 3 4 6 Pump diameter (mm) 57 44 38 44 Pump clearance level 1 2 2 1 Pump depth (m) 890 1377 1138 1430 Sucker-rod string (mm × m) 25∗890 19∗673+22∗704 19∗672+22∗466 19∗720+22∗710 middle depth of reservoir (m) 1000 1492 1569 1600 Crude oil density (kg/m3) 795 857 857 850 Water content (%) 95 97 92 98 Dynamic liquid level (m) 880 1347 757 1340 Casing pressure (Pa) 0.2 0.3 0.3 0.2 Oil pressure (Pa) 0.3 0.3 0.3 0.2 Fluid dynamic viscosity (Pa·s) 0.006 0.007 0.007 0.006 Gas oil ratio (m3/ m3) 20 19 40 80 Plunger length (m) 1.2 1.2 1.2 1.2 Clearance length (m) 0.5 0.5 0.5 0.5Figure 8
Dynamometer sensor.Figure 9
Suspension point dynamometer card comparison between the simulated and measured.
(a)
Well #1 (b)
Well #2 (c)
Well #3 (d)
Well #4


## 6. Dynamic Response Comparison

Pump dynamometer card is a closed graph recording plunger loads versus plunger displacement over a SRPS cycle. It is determined by multifactors such as stroke, stroke frequency, pump diameter, pump depth, dynamic liquid level, and the other oil well parameters. The difference between the current SRPS model and the improved SRPS model proposed in this paper is whether to consider the PFFP. In view of this, the oil operating status is divided into two forms as follows: (1) the fluid always keeps pace with the plunger when it is being sucked into the pump; (2) the fluid is unable to follow with the plunger. Then two oil wells are selected: well#1 with sufficient OWD and well#2 with insufficient OWD. Figure10 describes the plunger and fluid velocity during upstroke based on the improved SRPS model. It can be concluded that when the well has sufficient OWD and the fluid possesses good ability of keeping pace with the plunger; when the well has insufficient OWD, the fluid velocity lags behind the one of plunger at the beginning, whereas it is faster after a certain time period. Then the comparisons of pump dynamometer card are executed by the two models; meanwhile their individual pump fullness result is given at the end of the legend (Figure 11).Figure 10
Plunger and fluid velocity.
(a)
Well #1 (b)
Well #2Figure 11
Comparison of pump dynamometer card.
(a)
Well #1 (b)
Well #2Pump dynamometer card is very important method and normally used to diagnose the pump operations, particularly whose shape is more similar to a rectangle, the pump is closer to the fullness [28]. Hence, from Figure 11, it can be known that the pump fullness of current SRPS model is higher than that of the improved one depending on the qualitative judgment, and this conclusion is consistent with the quantitative calculation result, meanwhile this gap for the well #2 is bigger than well #1. According the above description, we can know that the pump fullness calculation result is on the high side if the PFFP is not considered, and this phenomenon will become more obvious for the well with insufficient OWD.The load presented by the left, upper, right, and lower borderline of pump dynamometer card is the results of pump moving from phase 1 to 4 in turn. Therefore, its upper borderline describes the pump load when fluid is pumped into the barrel. From Figure11(a), it can be found that the upper borderline load simulated by the improved SRPS model is significantly larger than the one of current. For Figure 11(b), this difference can be neglected. Based on (1) and (2), the pump pressure keeps constant asps when fluid is sucked into the pump and then the corresponding pump load is a fixed value. Due to the new model takes into account of PFFP, a pressure drop will be brought, which will vary with the fluid velocity. Combined with (1), this pressure drop will lead to the pump load increases. Figure 10 shows that the fluid velocity of well #1 is larger than that of well #2, and this is the reason why the difference of upper borderline load between the two models is obvious. This illustrates the error on the upper borderline load of pump dynamometer card depending on the velocity at which the fluid flows into the pump.


## 7. Optimization Test

Based on the above models, a simulation and optimization software is developed by MATLAB. One well with insufficient OWD is tested. Its original swabbing parameters are as follows: strokeS is 3 m, stroke frequencyns is 3 min−1, pump diameterDd is 57mm, pump depthLpd is 900 m, crank balance radiusrc is 1.2 m, and rod string combinationdj ×Lj is 25 mm × 524 m+22 mm × 376 m. And its swabbing parameters after optimizing are as follows strokeS is 3 m, stroke frequencyns is 2.5 min−1, pump diameterDd is 57mm, pump depthLpd is 990 m, crank balance radiusrc is 0.9 m, and rod string combinationdj ×Lj is 22 mm × 426 m+19 mm × 564 m. The comparisons before and after optimization are shown in Figure 12 and Table 2.Table 2
Specific parameters comparison before and after optimization.
Suspension point load Crank torque Motor input power Pump fullness Maximum Minimum Amplitude Standard deviation Balance degree Maximum Minimum mean Before optimization 55.5 kN 19.0 kN 36.5 kN 10.7 kN 85.4 % 15.9 kW 0.9 kW 6.2 kW 72.5 % After optimization 44.1 kN 13.2 kN 30.9 kN 9.4 kN 97.8 % 10.6 kW 1.0 kW 5.4 kW 90.8 % ↓ 20.5 % ↓ 30.5 % ↓ 15.3 % ↓ 12.1 % ↑ 14.5 % ↓ 33.3 % ↑ 11.1 % ↓ 12.9 % ↑ 25.2 %Figure 12
Dynamic response comparison before and after optimization.
(a)
Suspension point dynamometer card (b)
Crank torque (c)
Motor input powerFrom the above comparison results, some conclusions are obtained:(1) After optimization, the maximum and minimum of suspension point load are all decreased, and the load amplitude is lowered by 15.3 %. It can contribute to enhance the rod string life and prolong the maintenance period.(2) After optimization, the standard deviation of crank torque is reduced by 12.1 %, and its balance degree is raised by 14.5 %. It illustrates that the load torque fluctuation is cut down, which can minimize damage to transmission parts and improve the motor efficiency.(3) After optimization, the pump fullness is improved by 25.2 %. It is conducive to improve pump efficiency and production.(4) After optimization, it plays a role of peak shaving and valley filling for motor input power and extends the operational life, meanwhile power saving rate is 12.9 %.


## 8. Conclusions

(1) In this article, an improved SRPS model is presented, considering the couple effect of pumping fluid and plunger motion on the dynamic response of SRPS, instead of the existing models that assume the pumping fluid volume is always equal to plunger travelling volume. The Runge-Kutta method is applied to solve the whole system model trough transforming the rod string longitudinal vibration equation into ordinary differential equations. And the SRPS model’s precision has been validated by adopting surface dynamometer card.(2) Two oil wells are served to compare the difference between the current SRPS model and the improved one. The results indicate that the current SRPS model is relatively low in calculating pump fullness, and this gap will increase with the reduction of OWD; the influence of fluid flowing into the pump on pump load cannot be ignored when the fluid velocity is high. Therefore, the PFFP is necessary to be considered in order to improve the simulation accuracy.(3) On the basis of the improved SRPS model, a multitarget optimization model is proposed in purpose of improving production, decreasing load fluctuation and saving energy. By comparison, the optimal scheme can achieve the decreasing of maximum and minimum suspension point load, crank torque fluctuation and energy consumption, as well as improving the system balance degree and pump fullness. In summary, it improves the dynamic behavior of SRPS.



---

*Source: 4746210-2018-10-30.xml*