Exercise 8.2
Bo vs Pressure - Bubble-Point Peak
The oil formation volume factor Bo tells you how many reservoir barrels each stock-tank barrel occupies downhole, the number that turns OOIP into a tank gauge. For a saturated, black-oil reservoir Bo does something that trips up every junior engineer the first time they see it: it does not keep rising as pressure climbs. It peaks at the bubble point Pb and then bends back down.
The physics is a tale of two regimes for an OML 58 crude (API = 35, gamma_g = 0.70, T_F = 200 °F, Rs_at_Pb = 650 scf/STB):
- Below
Pb(saturated): drop the pressure and gas comes out of solution.
Less dissolved gas means a smaller, denser oil; Bo falls. Model the dissolved gas as a linear ramp Rs(P) = Rs_at_Pb * (P / Pb) and feed it to Standing: Bo = standing_Bo(Rs(P), gamma_g, gamma_o, T_F).
- Above
Pb(undersaturated): all the gas is already dissolved, so adding
pressure just compresses the liquid, so Bo shrinks slowly with oil compressibility c_o = 1.0e-5 /psi: Bo = Bo_pb exp(-c_o (P - Pb)), where Bo_pb is Bo evaluated right at Pb. The result: Bo is maximized at Pb.
Use the embedded api_to_sg, standing_bubble_point, and standing_Bo correlations (already in the starter). Then build:
bo_curve(P_array, API, gamma_g, T_F, Rs_at_Pb): return a NumPy array of
Bo for each pressure in P_array, applying the saturated branch where P <= Pb and the undersaturated branch where P > Pb. Compute gamma_o = api_to_sg(API) and Pb = standing_bubble_point(...) inside.
- A matplotlib figure of
Bo(y) vsP(x) using the supplied pressure grid
p_psia (which deliberately includes Pb itself), with axis labels, a title, and a vertical axvline at Pb marking the bubble point.
> Think about it: a reservoir produced from above Pb shows Bo rising > as you draw it down, right up until you cross the bubble point and free gas > appears. Why does that single inflection dominate material-balance and > surface-facility design?
Stuck? Reveal hints one at a time — they progress from nudge to near-solution.
visibilityReveal reference solutionexpand_more
Try solving it yourself first — the hints walk you through it. The solution below is one valid approach; yours may differ and still be correct.
import numpy as np
import matplotlib.pyplot as plt
# ── Embedded PVT correlations (Standing) - do not modify ────────────────────
def api_to_sg(API):
return 141.5 / (API + 131.5)
def standing_bubble_point(Rs, gamma_g, T_F, API):
exponent = 0.00091 * T_F - 0.0125 * API
return 18.2 * ((Rs / gamma_g) ** 0.83 * 10 ** exponent - 1.4)
def standing_Bo(Rs, gamma_g, gamma_o, T_F):
F = Rs * np.sqrt(gamma_g / gamma_o) + 1.25 * T_F
return 0.972 + 1.47e-4 * F ** 1.175
# ── OML 58 crude ────────────────────────────────────────────────────────────
API = 35.0
gamma_g = 0.70
T_F = 200.0
Rs_at_Pb = 650.0 # scf/STB dissolved at the bubble point
C_O_PER_PSI = 1.0e-5 # undersaturated oil compressibility
# Bubble point for this fluid - used to anchor the pressure grid.
PB_PSIA = standing_bubble_point(Rs_at_Pb, gamma_g, T_F, API)
# Pressure grid (psia). Built so that Pb is one of the grid points, so the
# curve's maximum lands exactly on the bubble-point pressure.
_below = np.linspace(200.0, PB_PSIA, 40)
_above = np.linspace(PB_PSIA, 4500.0, 40)[1:]
p_psia = np.concatenate([_below, _above])
def bo_curve(P_array, API, gamma_g, T_F, Rs_at_Pb):
"""Two-branch Bo(P): saturated below Pb, undersaturated above Pb.
Returns a NumPy array of Bo (RB/STB) for each pressure in P_array.
"""
P_array = np.asarray(P_array, dtype=float)
gamma_o = api_to_sg(API)
Pb = standing_bubble_point(Rs_at_Pb, gamma_g, T_F, API)
Bo_pb = standing_Bo(Rs_at_Pb, gamma_g, gamma_o, T_F)
# Saturated: dissolved gas ramps linearly with pressure up to Pb.
rs_p = Rs_at_Pb * (P_array / Pb)
bo_sat = standing_Bo(rs_p, gamma_g, gamma_o, T_F)
# Undersaturated: liquid simply compresses above Pb.
bo_undersat = Bo_pb * np.exp(-C_O_PER_PSI * (P_array - Pb))
return np.where(P_array <= Pb, bo_sat, bo_undersat)
bo_rb_stb = bo_curve(p_psia, API, gamma_g, T_F, Rs_at_Pb)
fig, ax = plt.subplots(figsize=(9, 5))
ax.plot(p_psia, bo_rb_stb, color="tab:green", lw=2, label="Bo(P)")
ax.axvline(PB_PSIA, color="crimson", ls="--", lw=1.5,
label=f"Pb = {PB_PSIA:.0f} psia")
ax.set_xlabel("Pressure (psia)")
ax.set_ylabel("Bo (RB/STB)")
ax.set_title("OML 58 oil FVF vs pressure - peak at the bubble point")
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
print("Pb (psia):", PB_PSIA)
print("Bo array:", bo_rb_stb)
lockCopying code is a Full Access feature.