libbsslo.bas

本系统是给大庆油田做的一个示例程序

Attribute VB_Name = "LibBSsLo"
Option Explicit
Sub Lopt_Cal(pb_well As Single, p As Single, twell As Single, nw As Single, sp As Single, SCrudeOil As Single, _
             lopt As Single)

    Dim e_sp As Single
    Dim so As Single, sw As Single, sng As Single
    Dim p0 As Single, T0 As Single
    Dim bwpt As Single
   
    Dim aa1 As Single, aa2 As Single
    Dim zpst As Single, sspst As Single, bopst As Single, lopst As Single
    Dim bopt As Single
    
    e_sp = 0
    Call BoSsLo(p, twell, nw, sp, SCrudeOil, pb_well, zpst, sspst, bopst, lopst)

    so = SCrudeOil
    sw = 1000
    sng = 1.0209
    
    p0 = 0.1013 * 1000000#
    T0 = 273.15

    bwpt = 1
    bopt = bopst

    aa1 = (1 - nw) * so + nw * sw + (1 - nw) * sp * sng
    
    If nw <= 0.97 Then
        aa2 = (1 - nw) * bopt + nw * bwpt + (1 - nw) * (1 - e_sp) * (sp - sspst) * p0 / T0 * twell * zpst / p
    Else
        aa2 = (1 - 0.97) * bopt + nw * bwpt + (1 - 0.97) * (1 - e_sp) * (sp - sspst) * p0 / T0 * twell * zpst / p
    End If
    lopt = aa1 / aa2
End Sub

Sub BoSsLo(p As Single, twell As Single, nw As Single, sp As Single, SCrudeOil As Single, pb_well As Single, _
               zpt As Single, sspt As Single, bopt As Single, lopt As Single)
    Dim e_sp As Single
    Dim so As Single, sw As Single, sng As Single
    Dim p0 As Single, T0 As Single
    Dim bwpt As Single, Bo1 As Single, Bo2 As Single, Bo As Single
    Dim Cb1 As Single, Cb2 As Single
    Dim aa1 As Single, aa2 As Single
    Dim alfas As Single '溶解度系数
    Dim ngas As Single
  
    e_sp = 0
    
    p0 = 0.1013 * 1000000#
    T0 = 273.15
    
    alfas = 3.5 * 0.000001
    If p <= pb_well Then
        sspt = alfas * (p - p0)
    Else
        sspt = alfas * (pb_well - p0)
    End If
    
    If sspt >= 0.9 * sp Then sspt = 0.9 * sp
           
    If p <= pb_well Then
        Bo1 = 1.022
        Cb1 = 8.2243 * 0.000000001
        Bo = Bo1 + Cb1 * p
    Else
        Bo2 = 1.11
        Cb1 = 2.5143 * 0.000000001
        Bo = Bo2 - Cb2 * (p - pb_well)
    End If
    If Bo <= 1 Then Bo = 1
    
    bopt = Bo
    bwpt = 1
    
    zpt = 0.9214
    Call Zcal(twell, p, ngas, zpt)
    
    so = SCrudeOil
    sw = 1000
    sng = 1.0209
    
    aa1 = (1 - nw) * so + nw * sw + (1 - nw) * sp * sng
    
    If nw < 0.97 Then
       aa2 = (1 - nw) * Bo + nw * bwpt + (1 - nw) * (1 - e_sp) * (sp - sspt) * p0 / T0 * twell * zpt / p
    Else
       aa2 = (1 - 0.97) * Bo + nw * bwpt + (1 - 0.97) * (1 - e_sp) * (sp - sspt) * p0 / T0 * twell * zpt / p
    End If
    lopt = aa1 / aa2
End Sub
Sub Zcal(T As Single, p As Single, ngas As Single, zgas As Single)
    Dim Xstart As Single, Xend As Single, Deta As Single, Eps As Single, Rootx As Single
    
    Xstart = 0.5
    Xend = 1.2
    Deta = 0.05
    Eps = 0.00001
    Call RootZ(Xstart, Xend, Deta, Eps, Rootx, T, p, ngas)
    zgas = Rootx
End Sub

Sub RootZ(Xstart As Single, Xend As Single, Deta As Single, Eps As Single, _
          Rootx As Single, T As Single, p As Single, ngas As Single)
    Dim XX As Single, X1 As Single, fx As Single, F1 As Single, X2 As Single
    Dim xb As Single, F2 As Single, X0 As Single, F0 As Single
    
    XX = 0
    X1 = Xstart
    XX = X1
    Call FuncZ(XX, fx, T, p, ngas)
    F1 = fx
    Do
       X2 = X1 + Deta
       xb = X1 + Deta
       XX = X2
       Call FuncZ(XX, fx, T, p, ngas)
       F2 = fx
       If F1 * F2 < 0 Then
          Exit Do
       Else
          X1 = xb
          F1 = F2
       End If
    Loop
    Do
       X0 = (X1 + X2) / 2
       XX = X0
       Call FuncZ(XX, fx, T, p, ngas)
       F0 = fx
       F2 = F0
       If Abs(X0 - X2) <= Eps Or Abs(F0) <= Eps Then
          Rootx = X0
          Exit Do
       ElseIf F1 * F2 > 0 Then
          F1 = F0
          X1 = X0
       ElseIf F1 * F2 < 0 Then
          X2 = X0
       End If
    Loop
End Sub
Sub FuncZ(XX As Single, fx As Single, T As Single, p As Single, ngas As Single)
    Dim tc As Single, pc As Single, tr As Single, pr As Single, c6 As Single, c7 As Single
  
    If ngas = 1 Then
       tc = 233.9863
       pc = 4.5022 * 1000000#
    Else
       tc = 132.5
       pc = 3.766 * 1000000#
    End If
    tr = T / tc
    pr = p / pc
    c6 = (0.3051 - 1.0467 / tr - 0.5783 / tr ^ 3) * 0.27 * pr / tr
    c7 = (0.5353 - 0.6123 / tr + 0.6818 / tr ^ 3) * (0.27 * pr / tr) ^ 2
    fx = XX ^ 3 - XX ^ 2 - c6 * XX - c7
End Sub
Rem 油气水三相混合液体积系数
Sub Bogw_Cal(pb_well As Single, ps As Single, p As Single, twell As Single, nw As Single, sp As Single, _
             SCrudeOil As Single, Bogw)
    Rem Bogw——油气水混合液的体积系数
   
    Dim p0 As Single, T0 As Single
    Dim zpst As Single, sspst As Single, bopst As Single, lopst As Single
    Dim zpt As Single, sspt As Single, bopt As Single, lopt As Single
    Dim bwpt As Single
    Dim x As Single
    
    p0 = 0.1013 * 1000000#
    T0 = 273.15
    
    Call BoSsLo(ps, twell, nw, sp, SCrudeOil, pb_well, zpst, sspst, bopst, lopst)
    Call BoSsLo(p, twell, nw, sp, SCrudeOil, pb_well, zpt, sspt, bopt, lopt)
    
    bwpt = 1
    x = (1 - nw) * bopt + nw * bwpt + (1 - nw) * (sp - sspst) * p0 / T0 * twell * zpt / p
    Bogw = 1 / x
    If Bogw <= 0.82 Then Bogw = 0.82
End Sub