home *** CD-ROM | disk | FTP | other *** search
/ Hall of Fame / HallofFameCDROM.cdr / oilfield / gasper.lzh / GASPER.BAK < prev    next >
Text File  |  1980-01-02  |  8KB  |  252 lines
  1. c*****************************************************************
  2. c          This program simulates gas well performance
  3. c    It computes either static of flowing bottom
  4. c    hole pressure for both sweet and/or sour gas
  5. c    wells 
  6.     character yes*5,no*5,ans*5,what*5,fname*15
  7.     real h(6),w(6),mu,pt
  8.     common theta,tem,w,pt,co2,h2s
  9.      &  ,h,q,tt,beta,e,dia,mu,sg,dept,ffact,zf
  10. c    read data for Gaussian six points quadrature
  11.     write(*,*)'PLEASE GIVE YOUR NAME IN QUOTES'
  12.         write(*,*)'LIKE '' YOUR NAME'', AND HIT ''ENTER'''
  13.     read(*,*)fname
  14.     write(*,*)' HENRY WELCOMES YOU,  ',fname
  15.     write(*,*)'HAVE YOU LOOK AT THE README FILE'
  16.     write(*,*)'********************************'
  17.     write(*,*)' JUST ANSWER YES OR NO THEN HIT ''ENTER'''
  18.     write(*,*)'********************************'
  19.     read(*,*)what
  20.     if(what.eq.'no')then
  21.     write(*,*)' PLEASE DO BEFORE RUNNING THIS PROGRAM'
  22.     write(*,*)'********************************'
  23.     write(*,*)'TO READ, AT THE PROMPT,ENTER '' TYPE README'''
  24.     stop
  25.     endif
  26.     write(*,1)
  27. 1    format(4x,'THIS PROGRAM COMPUTES STATIC OR FLOWING',//,
  28.      &  4x,'BOTTOM HOLE PRESSURES IN SWEET AND/OR SOUR',//,
  29.      &  4x,'GAS WELLS, TO USE, JUST ANSWER THE QUESTIONS',//
  30.      &  ,4x,'PLEASE DIRECT ALL INQUIRIES TO',//,4x,
  31.      &  ' HENRY A. OHEN , SCHOOL OF PETROLEUM',
  32.      &  //,4x,'& GEOLOGICAL ENGINEERING, UNIVERSITY',//,4x
  33.      &  , 'OF OKLAHOMA, NORMAN OK 73019',//)
  34.     open(10,file='gaus.dat')
  35.     do 2 i=1,6
  36.     read(10,*)h(i),w(i)
  37.  2    continue
  38.     imax=100
  39.     write(*,3)
  40.  3    format(4x,'IS YOUR GAS SWEET',//,4x,
  41.      &  ' ANSWER YES OR NO AND HIT ''ENTER''',//) 
  42.     read(*,*)ans
  43.     if(ans.eq.'no')then
  44.     write(*,4)
  45.  4    format(4x,' YOU MUST KNOW THE MOLE PERCENT OF CARBON DIOXIDE',
  46.      &  //,4x,' AND THE MOLE PERCENT OF HYDROGEN SULPHIDE',//,
  47.      & 'DO YOU HAVE THIS INFORMATION? YES OR NO AND HIT ''ENTER''',//)
  48.     read(*,*)ans
  49.     if(ans.eq.'no')then
  50.     write(*,5)
  51.  5    format(4x,' YOU ARE NOT READY,GET THE REQUIRED INFORMATION'
  52.      &  ,//)
  53.     stop
  54.     else
  55.     write(*,6)
  56.  6    format(4x,'GIVE MOLE FRACTION OF CO2 AND H2S IN THAT ORDER',/
  57.      & ,4x,' HIT ''ENTER'' FOR EACH ENRTY')
  58.     read(*,*)co2,h2s
  59.     endif
  60.     else
  61.     co2=0
  62.     h2s=0
  63.     endif
  64.     write(*,7)
  65.     read(*,*)dept
  66.     write(*,9)
  67.     read(*,*)tt
  68.     write(*,10)
  69.     read(*,*)pt
  70.     write(*,12)
  71.     read(*,*)dia
  72.     write(*,14)
  73.     read(*,*)beta
  74.     write(*,16)
  75.     read(*,*)theta
  76.     write(*,17)
  77.     read(*,*)sg
  78.     write(*,19)
  79.     read(*,*)tol
  80.     write(*,20)
  81.     read(*,*)q
  82.  7    format(4x,'GIVE WELL DEPTH(ft), HIT ''ENTER''')
  83.  9    format(4x,'TUBING HEAD TEMP.(of), HIT ''ENTER''')
  84.  10     format(4x,'GIVE PRESURE (psi), HIT ''ENTER''')
  85.  12     format(4x,'PIPE DIAMETER (inches), HIT ENTER''')
  86.  14     format(4x,'GIVE TEMPERATURE GRADIENT (of/ft), HIT ''ENTER''')
  87.  16     format(4x,'GIVE HORIZONTAL INCLINATION (degree), HIT ''ENTER''')
  88.  17     format(4x,'GIVE GAS GRAVITY HIT ''ENTER''')
  89.  19     format(4x,'GIVE TOLERANCE,1.0e-03,IS FINE, HIT ''ENTER''') 
  90.  20     format(4x, 'GIVE GAS FLOW RATE (MMscfd) , HIT ''ENTER''')
  91.  
  92.     theta=theta/180.0
  93.     e=0.0006
  94.     xguess=pt+0.3*pt
  95.  15    call fcn(xguess,f,fp)
  96.     xx=xguess-f/fp
  97.     if(ij.gt.imax)then
  98.     write(*,*)'ITERATION EXCEEDED IN MAIN LOOP, INCREASE TOL'
  99.     stop
  100.     endif
  101.     ch=(xx-xguess)
  102.     if(abs(ch).gt.tol)then
  103.     xguess=xx
  104.     ij=ij+1
  105.     go to 15
  106.     endif
  107.     pw=xx
  108.     open (13,file='res')
  109.     if(q.eq.0)then
  110.     write(13,*)'    ************************'
  111.     write(13,*)'    THIS IS A STATIC WELL'
  112.     write(13,*)'    ************************'
  113.     else
  114.     write(13,*)'    ************************'
  115.     write(13,*)'    THIS IS A FLOWING WELL'
  116.     write(13,*)'    ************************'
  117.     endif
  118.     write(13,31)zf,mu,ffact,pw
  119.  31    format(4x,'Z-FACTOR IS',4x,f7.4,//,
  120.      &  4x,'GAS VISCOSITY AT BOTTOM IS',4x,f7.4,//,
  121.      &   //,4x,'MOODY FRICTION FACTOR IS'
  122.      &  ,4x,f7.4,//4x,
  123.      &  'PRESURE AT GIVEN DEPT IS',2x,f13.5,//)
  124.     write(*,8)
  125.  8    format(//,2x,'  ****************************************')
  126.     write(*,*)'     FIND YOUR RESULTS IN A FILE CALLED ''res'''
  127.     write(*,*)'     ****************************************'
  128.     write(*,*)'     TO READ ENTER ''TYPE RES'''
  129.     write(*,33)fname
  130.     write(*,*)'****************************'
  131. 33    format(//,4x,'HAVE A GOOD DAY,  ',a15,)
  132.     end
  133. c******************************************************
  134.     subroutine fcn(x,f,fp)
  135.     real h(6),w(6),x,mu,integral
  136.     common theta,tem,w,pt,co2,h2s,
  137.      &  h,q,tt,beta,e,dia,mu,sg,dept,ffact,zf
  138.     tem=tt+beta*dept+460
  139.     call solve(x,integral,fp)
  140.     f=1000*dept-integral
  141.     return
  142.     end
  143.     
  144. c***************************************************
  145.     subroutine solve(p,res,fp)
  146.     real t,res,mu,dia,sg,h(6),w(6)
  147.     common theta,tem,w,pt,co2,h2s
  148.      &  ,h,q,tt,beta,e,dia,mu,sg,dept,ffact,zf
  149.     shc=co2*h2s
  150.     rl=dept/(sin(theta))
  151.     ee=120*(shc**0.9-shc**1.6)+15*(h2s**0.5-h2s**4.0)
  152.     tc=170.491+307.344*sg-ee
  153.     pc=709.604-58.715*sg
  154.     pc=pc*tc/((tc+ee)+h2s*(1.0-h2s)*ee)
  155.     ttc=tem/tc
  156.     pr=p/pc
  157.     call zfact(ttc,pr,z)
  158.     zf=z
  159.     call gvis(ttc,tc,pr,z,pc,sg,fm)
  160.     mu=fm
  161.     if(q.ne.0)then
  162.     ren=20000*q*sg/mu/dia
  163.     call fact(ren,fr)
  164.     ffact=fr
  165.     endif
  166.     a1=(53.34*p/sg/tem/zf+111.1*q*q/(dia**4)/p)
  167.     a2=(2.6665*fr*q*q/4)/(dia**5)
  168.     a3=p*p/tem/tem/1000/zf/zf
  169.     fp=-a1/(a2+a3)
  170.     call gaus(fr,z,p,res)
  171.     return
  172.     end
  173.  
  174.     subroutine fact(ren,f)
  175. c******************************************************
  176. c     Subroutine calculates the friction factor by solving
  177. c     The Colebrooks and White Implicit Equation for friction factor
  178. c     using the Newton Raphson's method
  179. c*****************************************************
  180.     real f,ren,b,fn,ff,h(6),w(6)
  181.     common theta,tem,w,pt,co2,h2s
  182.      &  ,h,q,tt,beta,e,dia,mu,sg,dept,ffact,zf
  183. c     Assume Laminar Flow as first approximation
  184.     fo=64/ren
  185.     do 2 i=1,50
  186. c     Start Newton Iteration
  187.     ff=1.74-0.8686*alog(2*e/dia+18.7/sqrt(fo)/ren)-1.0/sqrt(fo)
  188.     aa=(2*e/dia+18.7/ren/sqrt(fo))
  189.     fn=0.5*(1.0+16.24283/aa/ren)/(fo**1.5)
  190.         f=fo-ff/fn
  191.     if(abs((f-fo)/fo).le.1e-05)go to 5
  192.     fo=f
  193.  2    continue
  194.     write(*,*)'iteration failed in ffact after ',i,'times'
  195.     stop
  196.  5    return
  197.     end
  198. c************************************************************
  199.     subroutine gaus(fr,z,xy,res)
  200. c************************************************************
  201.     real y(6),h(6),w(6),length
  202.     common theta,tem,w,pt,co2,h2s
  203.      &  , h,q,tt,beta,e,dia,mu,sg,dept,ffact,zf
  204.     a=pt
  205.     b=xy
  206.     do 20 i=1,6
  207.     y(i)=(h(i)+1.0)*(b-a)/2.0+a
  208.  20    continue
  209.     sum=0
  210.     length=dept/sin(theta)
  211.     do 30 i=1,6
  212.     a1=y(i)*53.34/sg/tem/z
  213.     a2=111.1*q*q/(dia**4)/y(i)
  214.     a3=2.6665*fr*q*q/4/(dia**5)
  215.     a4=y(i)*y(i)/tem/z/z/tem/1000
  216.     g=(a1+a2)/(a3+a4)
  217.     sum=sum+w(i)*g*0.5*(b-a)
  218.  30    continue
  219.     res=sum
  220.     return
  221.     end
  222.  
  223.         subroutine zfact(tr,pr,z)
  224. c       this routine calculates z-factor using  Standing modification
  225. c       to the Brill and Beggs correlation for curve-fitting the Stand-
  226. c       ing and Katz reduced pressure-reduced temperature z-factor chart
  227. c
  228.         a=1.39*(tr-0.92)**0.5-0.36*tr-0.101
  229.         b=(0.62-0.23*tr)*pr
  230.         c=(0.066/(tr-0.86)-0.037)*pr**2
  231.         d=(0.32/(10.**(9.*(tr-1.))))*pr**6
  232.         e=b+c+d
  233.         f=(0.132-0.32*alog10(tr))
  234.         g=10.**(0.3106-0.49*tr+0.1824*tr**2)
  235. c
  236.         z=a+(1.-a)*exp(-e)+f*pr**g
  237.         return
  238.         end
  239. c**********************************************************
  240.  
  241.     subroutine gvis(tr,tc,pr,z,pc,sg,fm)
  242. c     Subroutine computes gas viscosity using the 
  243. c    modified Lee et al correlation
  244.     dg=(2.7*pr*pc*sg)/(z*tr*tc)/(16.10846e03)
  245.     ak=((9.4+0.58*sg)*(tr*pc)**(1.5))/(209.0+
  246.      &  551*sg+tr*tc)
  247.     xx=3.5+986/tc/tr+0.29*sg
  248.     yy=2.4-0.2*xx
  249.     fm=ak*exp(xx*dg**yy)*1e-04
  250.     return
  251.     end
  252.