Exercise 19.3
The Over-Booking Trap - Capped vs Uncapped EUR
Project 2's wells decline steeply enough that the terminal-decline cap never bites: confirm that removing it barely moves the field EUR. Then build the trap deliberately: take a long-life well (qi = 1500, Di = 0.008/month) and sweep the Arps exponent b from 0.3 toward 1.0, computing EUR with the cap (modified_hyperbolic) and without it (plain hyperbolic). Show how the uncapped EUR explodes as b → 1 while the capped EUR stays bounded, and explain, in reserves-audit terms, why booking the uncapped number is indefensible: what does a b near 1 imply about a tail that is still flowing a century out?
---
The terminal-decline cap in Project 2 is not decoration; it is the single most important reserves safeguard in the chapter. This exercise shows you why, by removing it and watching what happens.
There are two parts. First, confirm the honest finding: Project 2's real wells (FITS) decline steeply enough that removing the cap barely changes the field EUR. Then build the trap on purpose: a long-life well (TRAP_QI = 1500, TRAP_DI = 0.008/month) whose Arps exponent b you sweep toward 1.0, and watch the uncapped EUR run away while the capped EUR stays bounded.
The data, decline machinery, FITS, and constants (ECON, DAYS, DMIN, HORIZON, TRAP_QI, TRAP_DI, B_LADDER) are in the do-not-edit block.
Write two functions and two module-level values:
def eur(qi, Di, b, cap=True):
'''EUR in barrels: integrate the rate above the economic limit.
cap=True -> modified_hyperbolic (terminal switch); cap=False -> raw hyperbolic.'''
def field_inflation(fits):
'''Percent EUR inflation from removing the cap, for a set of wells {id:(qi,Di,b)}.'''Exact procedure:
- In
eur, taket = np.arange(0, HORIZON), buildqwith
modified_hyperbolic(t, qi, Di, b, DMIN) when cap else hyperbolic(t, qi, Di, b), keep q = q[q >= ECON], and return float(np.sum(q * DAYS)).
field_inflation(fits)returns(uncapped_total / capped_total - 1) * 100, where
each total sums eur(*fits[w], cap=...) over the wells. Then field_inflation_pct = field_inflation(FITS), expect ≈ 0 for the steep real wells, but the function must work for any fits (a long-life high-b field is not safe).
trap= a dict{b: {'capped': ..., 'uncapped': ..., 'inflation_pct': ...}}
for each b in B_LADDER, with inflation_pct = (uncapped / capped - 1) * 100.
Expose exactly these names: eur, field_inflation, field_inflation_pct, trap.
> Think about it: at b = 0.99 the uncapped forecast books roughly three > times the capped EUR: barrels the well would only "produce" by flowing above > 20 bopd for centuries. A reserves auditor caps the decline for exactly this > reason. What does a fitted b pinned at its upper bound tell you about whether > the hyperbolic model itself is the right one for that well?
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 io
import numpy as np
import pandas as pd
from scipy.optimize import curve_fit
import warnings
warnings.simplefilter("ignore")
# ── Verified Chapter 19 production data + decline machinery (do not edit) ─
def hyperbolic(t, qi, Di, b):
return qi / np.power(1 + b * Di * t, 1.0 / b)
def modified_hyperbolic(t, qi, Di, b, Dmin):
"""Arps hyperbolic that switches to a terminal exponential decline at Dmin."""
t_sw = (Di / Dmin - 1) / (b * Di) if (b > 0 and Di > Dmin) else 1e9
return np.where(t <= t_sw, hyperbolic(t, qi, Di, b),
hyperbolic(t_sw, qi, Di, b) * np.exp(-Dmin * (t - t_sw)))
PROD = """well,date,oil_bopd
OD-001,2025-01,2347.9
OD-001,2025-02,2278.1
OD-001,2025-03,2276.2
OD-001,2025-04,2200.2
OD-001,2025-05,2151.0
OD-001,2025-06,2000.3
OD-001,2025-07,1827.7
OD-001,2025-08,1807.2
OD-001,2025-09,1657.7
OD-001,2025-10,1674.6
OD-001,2025-11,1570.1
OD-001,2025-12,1561.9
OD-001,2026-01,1532.1
OD-001,2026-02,1510.5
OD-001,2026-03,
OD-001,2026-04,1280.1
OD-001,2026-05,1333.5
OD-001,2026-06,1158.4
OD-001,2026-07,1224.0
OD-001,2026-08,1103.0
OD-001,2026-09,1138.6
OD-001,2026-10,1021.0
OD-001,2026-11,1095.7
OD-001,2026-12,965.1
OD-003,2025-01,3081.7
OD-003,2025-02,3007.0
OD-003,2025-03,2792.2
OD-003,2025-04,-200.0
OD-003,2025-05,2431.6
OD-003,2025-06,2247.2
OD-003,2025-07,2180.3
OD-003,2025-08,2171.0
OD-003,2025-09,1979.1
OD-003,2025-10,1788.0
OD-003,2025-11,1694.2
OD-003,2025-12,1576.2
OD-003,2026-01,1568.3
OD-003,2026-02,1428.9
OD-003,2026-03,1270.9
OD-003,2026-04,1315.6
OD-003,2026-05,1207.4
OD-003,2026-06,1177.8
OD-003,2026-07,1074.7
OD-003,2026-08,960.8
OD-003,2026-09,968.1
OD-003,2026-10,899.4
OD-003,2026-11,806.9
OD-003,2026-12,770.1
OD-005,2025-01,1807.2
OD-005,2025-02,1737.1
OD-005,2025-03,1720.5
OD-005,2025-04,1644.8
OD-005,2025-05,1593.0
OD-005,2025-06,1539.4
OD-005,2025-07,1535.7
OD-005,2025-08,15000.0
OD-005,2025-09,1450.4
OD-005,2025-10,1329.5
OD-005,2025-11,1324.6
OD-005,2025-12,1331.2
OD-005,2026-01,1266.3
OD-005,2026-02,1234.4
OD-005,2026-03,1129.9
OD-005,2026-04,1123.2
OD-005,2026-05,1135.1
OD-005,2026-06,1080.7
OD-005,2026-07,1084.6
OD-005,2026-08,1061.8
OD-005,2026-09,1038.6
OD-005,2026-10,1004.2
OD-005,2026-11,885.3
OD-005,2026-12,890.1
OD-007,2025-01,2871.6
OD-007,2025-02,2862.5
OD-007,2025-03,2779.9
OD-007,2025-04,2584.2
OD-007,2025-05,2406.8
OD-007,2025-06,2317.7
OD-007,2025-07,2126.1
OD-007,2025-08,2064.2
OD-007,2025-09,1929.2
OD-007,2025-10,1911.9
OD-007,2025-11,1777.9
OD-007,2025-12,1671.4
OD-007,2026-01,1627.4
OD-007,2026-02,1509.9
OD-007,2026-03,1353.1
OD-007,2026-04,1332.1
OD-007,2026-05,1272.4
OD-007,2026-06,1204.4
OD-007,2026-07,1175.7
OD-007,2026-08,1082.7
OD-007,2026-09,1110.0
OD-007,2026-10,1099.1
OD-007,2026-11,953.7
OD-007,2026-12,980.7
"""
prod = pd.read_csv(io.StringIO(PROD))
ECON, DAYS, DMIN = 20.0, 30.4, 0.06 / 12 # economic limit, days/month, 6%/yr terminal decline
HORIZON = 24000 # months -- long enough for an uncapped b->1 tail
def clean_mask(q):
med = np.nanmedian(q)
return np.isfinite(q) & (q > 0) & (q < 3 * med)
# Project 2's real wells, fit once. FITS[well] = (qi, Di, b).
FITS = {}
for _w in sorted(prod.well.unique()):
_q = prod[prod.well == _w].oil_bopd.values.astype(float)
_ok = clean_mask(_q)
_t = np.arange(len(_q))
_popt, _ = curve_fit(hyperbolic, _t[_ok], _q[_ok], p0=[_q[_ok][0], 0.05, 0.5],
bounds=([0, 0, 0.01], [8000, 1, 0.99]), maxfev=5000)
FITS[_w] = (float(_popt[0]), float(_popt[1]), float(_popt[2]))
# The "trap" type well and the b values to sweep.
TRAP_QI, TRAP_DI = 1500.0, 0.008
B_LADDER = (0.3, 0.6, 0.9, 0.99)
# ── end do-not-edit ──────────────────────────────────────────────────────
def eur(qi, Di, b, cap=True):
"""EUR in barrels: integrate the rate above the economic limit. With cap=True
use the terminal-decline switch; with cap=False use the raw hyperbolic tail."""
t = np.arange(0, HORIZON)
q = modified_hyperbolic(t, qi, Di, b, DMIN) if cap else hyperbolic(t, qi, Di, b)
q = q[q >= ECON]
return float(np.sum(q * DAYS))
# (1) Removing the cap barely moves a STEEP field but explodes a long-life one.
def field_inflation(fits):
"""Percent by which removing the terminal-decline cap inflates the total EUR
of a set of wells {id: (qi, Di, b)}. ~0 for steep wells (the cap never bites);
large for long-life, high-b wells (the uncapped tail runs away)."""
capped = sum(eur(*fits[w], cap=True) for w in fits)
uncapped = sum(eur(*fits[w], cap=False) for w in fits)
return (uncapped / capped - 1) * 100
field_inflation_pct = field_inflation(FITS) # ~0.0 -- Project 2's real wells are safe
# (2) The trap: as b -> 1 the uncapped EUR explodes while the capped one stays bounded.
trap = {}
for b in B_LADDER:
c = eur(TRAP_QI, TRAP_DI, b, cap=True)
u = eur(TRAP_QI, TRAP_DI, b, cap=False)
trap[b] = dict(capped=c, uncapped=u, inflation_pct=(u / c - 1) * 100)
print(f"real field: removing the cap inflates EUR by {field_inflation_pct:.1f}% (the wells are safe)")
for b, v in trap.items():
print(f" trap b={b}: capped {v['capped']/1e6:.2f} uncapped {v['uncapped']/1e6:.2f} "
f"inflation {v['inflation_pct']:.0f}%")
lockCopying code is a Full Access feature.