计算井筒多相流的BeggsBrill方法讲解.docx
《计算井筒多相流的BeggsBrill方法讲解.docx》由会员分享,可在线阅读,更多相关《计算井筒多相流的BeggsBrill方法讲解.docx(32页珍藏版)》请在冰豆网上搜索。
计算井筒多相流的BeggsBrill方法讲解
作业一:
计算井筒多相流的Beggs-Brill方法
基本参数:
套管尺寸:
51/2’’,油管:
21/2’’,GOR=24
/
=0.9659,
=100cp,井底温度63℃,井口温度25℃,井口压力1MPa,产液137方/天,油层深度1490m,井底流压8MPa,含水率29%,泵的单级扬程4.6m。
设计步骤:
截图:
SubBeggsBrill(calAsSingle,rwAsSingle,roAsSingle,rgAsSingle,yco2AsSingle,yh2sAsSingle,yn2AsSingle,ynaclAsSingle,_ppcAsSingle,tpcAsSingle,d2AsSingle,d1AsSingle,erasAsSingle,lonAsSingle,alfa0AsSingle,T1AsSingle,pcal1AsSingle,T2AsSingle,qwAsSingle,qoAsSingle,QscAsSingle,daltalAsSingle,PCAL2AsSingle)
Do
Iflsum+daltal>lonThen
daltal=lon-lsum
EndIf
lmid=lsum+daltal/2
T=T1+(T2-T1)*lmid/lon
Do
daltap0=daltap
p=p1+Sgn(cal)*daltap/2'p利用假设的daltap计算,带入以后的p
lsum=lsum+daltal
p2=p1+Sgn(cal)*daltap
p1=p2
LoopUntillsum>=lon
PCAL2=p2
EndSub
psc=0.101325
tsc=293
densc=1.204*rg
deno=1000*ro
denw=1000*rw
qow=qo+qw
Ifqo>0Then
gor=Qsc/qo
CallPVTO(ro,rg,gor,T,p,pbo,rso,Bo,uo,sto)
EndIf
Ifqw>0Then
gwr=Qsc/qw
CallPVTW(T,p,ynacl,gwr,rsw,bw,uw,stw)
EndIf
Ifqow>0Then
GLR=Qsc/qow
RS=(rso*qo+rsw*qw)/qow
IfRS>=GLRThen
RS=GLR
EndIf
qgsc=qow*(GLR-RS)
Else
qgsc=Qsc
EndIf
Ifqgsc>0Then
ppr=p/ppc
tpr=T/tpc
Calldprzf(ppr,tpr,dr,z,cpr)
bg=psc/tsc*z*T/p
Calllgeug(rg,yn2,yco2,yh2s,p,T,z,ug)
qg=qgsc*bg
wg=densc*qgsc
densg=wg/qg
Else
qg=0
densg=0
bg=0
EndIf
Ql=qo*Bo+qw*bw
wl=deno*qo+denw*qw+densc*RS*qow
IfQl>0Then
densl=wl/Ql
fo=qo*Bo/Ql
EndIf
ul=uo*fo+uw*(1-fo)
stal=sto*fo+stw*(1-fo)
EndSub
gor=gorin*5.615
T=tin*1.8-460
p=145*pin
api=141.5/ro-131.5
Ifapi<=30Then
a1=0.0362
b1=1.0937
c1=25.724
Else
a1=0.0178
b1=1.187
c1=23.931
EndIf
tp=T+460
pbo=(gor/(a1*rg*Exp(c1*api/tp)))^(1/b1)
rso=a1*rg*p^b1*Exp(c1*api/tp)'ft^3/bbl
Ifp>=pboOrrso>=gorThen
rso=gor
EndIf
Ifapi<=30Then
a1=1.761*10^(-5)
b1=0.0004677
c1=-1.1811*10^(-8)
Else
a1=0.000011
b1=0.000467
c1=1.377*10^(-9)
EndIf
tp=T-60
Bo=1+a1*tp*(api/rg)+(b1+c1*tp*(api/rg))*rso
co=(-1433+5*rso+17.2*T-1180*rg+12.61*api)/1000001/p
Ifp>=pboOrrso>=gorThenBo=Bo*Exp(co*(pbo-p))
zz=3.0324-0.0203*api
yy=10^zz
xx=yy*T^(-1.163)
uod=10^xx-1
aa=10.715*(rso+100)^(-0.515)
bb=5.44*(rso+150)^(-0.338)
'uo=aa*uod^bb
uo=100'
Ifpbo>0And(p>=pboOrrso>=gor)Then
uo=uo*(p/pbo)^(2.6*p^1.187*Exp(-0.0000898*p-11.513))
EndIf
sto=(42.4-0.047*T-0.267*api)*Exp(-0.0007*p)/1000
pbo=pbo/145
rso=rso/5.615
EndSub
T=tin*1.8-460'F
p=145*pin
a1=2.12+0.00345*T-0.0000359*T^2
b1=0.0107-0.0000526*T+1.48*10^(-7)*T^2
c1=-8.75*10^(-7)+3.9*10^(-9)*T-1.02*10^(-11)*T^2
scnacl=1-(0.0753-0.000173*T)*ynacl
rws=(a1+b1*p+c1*p^2)*scnacl/5.615
Ifrws>gwrThen
rws=gwr
EndIf
Ifrws<=0Then
a1=0.9947+0.0000058*T+1.02*10^(-6)*T^2
b1=-4.228*10^(-6)+1.8376*10^(-8)*T-6.77*10^(-11)*T^2
c1=1.3*10^(-10)-1.3855*10^(-12)*T+4.285*10^(-15)*T^2
Else
a1=0.9911+6.350001*10^(-5)*T+8.5*10^(-7)*T^2
b1=-1.093*10^(-6)-3.497*10^(-9)*T+4.75*10^(-12)*T^2
c1=-5*10^(-11)+6.429*10^(-13)*T-1.43*10^(-15)*T^2
EndIf
tp=T-60
tl=5.1*10^(-8)*p+(5.47*10^(-6)-1.95*10^(-10)*p)*tp
scfm=(tl+(-3.23*10^(-8)+8.5*10^(-13)*p)*tp^2)*ynacl+1
bw=(a1+b1*p+c1*p^2)*scfm
tl=(T^0.5-0.0135*T)*(0.00276*ynacl-0.000344*ynacl^1.5)
scvis=1-0.00187*ynacl^0.5+0.000218*ynacl^2.5+tl
scpr=1+3.5*10^(-12)*p^2*(T-40)
uw=0.02414*scvis*scpr*10^(446/(T+208))
ax=1-6.94829*10^(-4)*p+8.878355*10^(-7)*p^2
ax=ax+2.36093*10^(-13)*p^4-6.42967*10^(-10)*p^3
ax=ax-4.19405*10^(-17)*p^5+2.86787*10^(-21)*p^6
stw=ax*(0.1746-0.00021*T)/2.2046
EndSub
'程序功能:
计算当地条件下的有关参数的过程
ga=9.807
don=d2-d1'
aer=3.14159*(d2^2-d1^2)/4
Vsl=Ql/86400/aer'm/s’
vsg=qg/86400/aer'
Vm=Vsl+vsg'
na=Vsl/Vm'
nr=Vm^2/(ga*don)'
IfVsl>0Then
nv=Vsl*(densl/(ga*stal))^0.25'
EndIf
denn=na*densl+(1-na)*densg'
un=(na*ul+(1-na)*ug)/1000
ren=denn*Vm*don/un'
fmn=(1/(1.14-2*Log(eras/don+21.25/ren^0.9)/Log(10)))^2'
EndSub
Ifna>0Then
ls1=316*na^0.302
ls2=0.0009252*na^(-2.4684)
ls3=0.1*na^(-1.4516)
ls4=0.5*na^(-6.738)
Ifna<0.00001Then
ident=5
ElseIfna>0.99999Then'
ident=6
ElseIfna<0.01Andnr=0.01Andnrident=1
ElseIf(na>=0.01Andna<0.4)And(nr>ls3Andnr<=ls1)Orna>=0.4And(nr>ls3Andnr<=ls4)Then
ident=2
ElseIfna<0.4Andnr>=ls1Orna>=0.4Andnr>=ls4Then
ident=3
ElseIfna>=0.01And(nr>=ls2Andnr<=ls3)Then
ident=4
EndIf
Else
ident=5
EndIf
EndSub
Ifalfa>=0Then
SelectCaseident
Case1
dd=0.011
ee=-3.768
ff=3.539
gg=-1.614
Case2
dd=2.96
ee=0.305
ff=-0.4473
gg=0.0978
Case3
dd=1
ee=0
ff=0
gg=0
EndSelect
Else
dd=4.7
ee=-0.3692
ff=0.1244
gg=-0.5056
EndIf
cc=(1-na)*Log(dd*na^ee*nv^ff*nr^gg)'
Ifcc<0Then
cc=0
EndIf
tm=Sin(1.8*alfa)
hui=1+cc*(tm-tm^3/3)'
Ifhui<0Then
hui=0
EndIf
SelectCaseident
Case1
dd=0.98
ee=0.4846
cc=0.0868
Case2
dd=0.845
ee=0.5351
cc=0.0173
Case3
dd=1.065
ee=0.5824
cc=0.0609
EndSelect
hl0=dd*na^ee*nr^(-cc)'
Ifhl0hl0=na
EndIf
hlsta=0.918*hl0*hui
Ifhlsta<0Then
hlsta=0
EndIf
Ifhlsta>1Then
hlsta=1
EndIf
EndSub
ga=9.807
SelectCaseident
Case4
ls2=0.0009252*na^(-2.4684)
ls3=0.1*na^(-1.4516)
ident=1
CallHOLDUP(alfa,ident,na,nv,nr,hlsta)
hlsta1=hlsta
ident=2
CallHOLDUP(alfa,ident,na,nv,nr,hlsta)
hlsta2=hlsta
tm=(ls3-nr)/(ls3-ls2)
hlsta=tm*hlsta1+(1-tm)*hlsta2
Case5
hlsta=0
Case6
hlsta=1
CaseElse
CallHOLDUP(alfa,ident,na,nv,nr,hlsta)
EndSelect
Ifhlstahlsta=na
EndIf
denm=densl*hlsta+densg*(1-hlsta)
tm=denm*ga*Sin(alfa)
Ifhlsta=0Orhlsta=1Then
sb=0
Else
yb=na/hlsta^2
yr=Log(yb)
Ifyb>1Andyb<1.2Then
sb=Log(2.2*yb-1.2)
Else
sb=yr/(-0.0523+3.182*yr-0.8725*yr^2+0.01853*yr^4)'=======16
EndIf
EndIf
fm=fmn*Exp(sb)
Tf=fm*denn*Vm^2/(2*don)
coef=1-denm*Vm*vsg/(1000000*p)
T1=(tm+Tf)/(1000000*coef)
daltap=daltal*T1
EndSub
a1=0.31506237
a2=-1.0467099
a3=-0.57832729
a4=0.53530771
a5=-0.61232032
a6=-0.10488813
a7=0.681570001
a8=0.68446549
a=a1+a2/tr+a3/tr^3
B=a4+a5/tr
C=a5*a6*tr
D=a7/tr^3
e=a8
G=-0.27*pr/tr
dr=-G
Do
ed=e*dr^2
exa=Exp(-ed)
f=G/dr+1+a*dr+B*dr^2+C*dr^5
f=f+D*dr^2*(1+ed)*exa
fp=-G/dr^2+a+2*B*dr+5*C*dr^4
fp=fp+2*D*dr*(1+ed-ed^2)*exa
dr=dr-f/fp
Ifdr=0Ordr=1Then
ExitSub
EndIf
LoopUntilAbs(f/fp)<0.00001
z=0.27*pr/(dr*tr)
ed=e*dr^2
exa=Exp(-ed)
dzdr=a+2*B*dr+5*C*dr^4
dzdr=dzdr+2*D*dr*(1+ed-ed^2)*exa
cr=1/pr-0.27/(z^2*tr)*dzdr/(1+dr*dzdr/z)
EndSub
mg=28.97*rg
R=0.008314
g1=yh2s*(5.7*rg-1.7)*0.00001
g2=yco2*(5*rg+1.7)*0.00001
g3=yn2*(5*rg+4.7)*0.00001
tr=1.8*T
G=(9.399999+0.02*mg)*tr^1.5
G=0.0001*G/(209+19*mg+tr)+g1+g2+g3
j=1.7-197.2/tr-0.002*mg
h=(3.5+986/tr+0.01*mg)*0.001^j
den=mg*p/(z*R*T)
ug=G*Exp(h*den^j)
EndSub
作业二抽油杆柱设计
基本参数配产:
30m3/d,泵深:
2300,动液面:
950m,冲程:
3.6,冲次6次/min,泵径:
56mm,油管内径:
73mm,地层水粘度:
1,8mPa.s,地层水密度:
1050kg/m3,抽油杆密度:
7850kg/m3,最小杆径:
19mm,杆使用系数:
0.95
设计步骤:
截图:
DimhcAsSingle'沉没度,m
DimSFAsSingle'杆使用系数
DimsAsSingle'冲程,m
DimqlAsSingle'产液量,m^3/d
DimdpAsSingle'泵径,mm
DimdeAsSingle'柱塞与泵筒间隙平均值,m
DimuwAsSingle'地层水粘度,mPa.s
DimnAsSingle'冲次,次/min
DimrouwAsSingle'地层水密度,kg/m^3
DimrousAsSingle'钢的密度,kg/m^3
DimtAsSingle'抽油杆最小抗张强度,MPa
DimPlAsSingle,PLiAsSingle,deltaPLAsSingle'输入应力范围比、计算应力范围比、两级抽油杆顶端面应力范围比的最大允许差值
DimdtAsSingle'油管内径,mm
DimfpAsSingle,frAsSingle,foAsSingle'柱塞截面积,抽油杆截面积,阀孔截面积,m^2
DimnrodAsInteger'杆级数
DimL0AsSingle,L1AsSingle'计算所得杆柱长度,m已设计出的杆柱总长度,m
DimdfoAsSingle'阀孔径,m
Dimqr()AsSingle'每米抽油杆质量,kg/m
DimLi()AsSingle'每级杆长度,m
Dimdr()AsSingle'抽油杆直径,mm
DimbAsSingle'失重系数
DimmAsSingle'油管内径与抽油杆直径之比
DimvmaxAsSingle'抽油杆柱最大下行速度,m/s
DimvpAsSingle,vfAsSingle'柱塞最大运动速度、阀孔液流速度,m/s
DimNReAsSingle,ksiAsSingle'雷诺数,阀流量系数
DimPmax0AsSingle,Pmin0AsSingle,PL0AsSingle'抽油杆上下冲程底端载荷,N,底端应力