A Python library for post-well analysis of rotary drilling data: Mechanical Specific Energy (MSE), d-exponent, bit hydraulics, founder point detection per formation and a WOB road map whose ROP projection comes only from the data. Stick-slip is analyzed on time-indexed windows, not on the depth log.
All results in this repository come from a synthetic benchmark well generated by
src/data/generator.py(fixed seed, reproducible). No field data has been run through the pipeline yet. See Limitations.
The library works on a depth-indexed surface log (depth, ROP, WOB, RPM, torque, flow, SPP, gamma) and, for torsional dynamics, on short time windows of downhole RPM and surface torque. It provides:
- MSE (Teale, 1965) and the efficiency comparison against rock strength proposed by Dupriest and Koederitz (2005).
- d-exponent and dc-exponent (Jorden and Shirley, 1966; Rehm and McClendon, 1971) and an Eaton-type pore pressure estimate.
- Bit hydraulics: TFA, nozzle velocity, bit pressure drop, HSI and jet impact force.
- Founder point detection per formation: hinge (two-segment) regression of ROP on WOB, with a significance test and a bootstrap interval for the breakpoint.
- WOB road map per formation: a WOB setpoint inside the tested range and the top-drive torque limit, with a projected ROP taken from the fitted ROP-WOB response. The projection can be negative and is zero when WOB was not varied.
- Stick-slip analysis on time windows: Stick-Slip Severity Index (SSI) from downhole RPM, FFT and Welch spectra, and a quarter-wave estimate of the string's torsional frequency.
Teale (1965), in oilfield units:
with WOB in lbf,
Dupriest and Koederitz (2005) scale MSE by a mechanical efficiency factor so that, when drilling efficiently, it approximates rock strength:
A rise in MSE with no change in rock strength points to energy lost to dysfunctions (bit balling, vibrations, worn cutters).
In the founder analysis, MSE is always recomputed from raw channels, and an off-bottom torque estimate (the zero-WOB intercept of a torque-vs-WOB fit) is subtracted from surface torque. On a rig this value is measured off bottom at connections.
Dupriest and Koederitz define the founder as the WOB at which the ROP response stops being linear. Below it ROP rises in proportion to WOB at roughly constant MSE; above it, extra WOB buys little ROP and MSE climbs.
FounderPointAnalyzer analyzes one formation (or one depth window) at a time:
- Keeps on-bottom samples and, by default, scales ROP to the interval's median RPM (assuming ROP proportional to RPM).
- Fits a straight line and a continuous hinge model, searching the breakpoint
$W_b$ on a grid between the 10th and 90th WOB percentiles:$$\text{ROP} = c + s_1 \cdot \text{WOB} + (s_2 - s_1) \cdot \max(0, \text{WOB} - W_b)$$ - Reports a founder at
$W_b$ only if the hinge beats the line with$F > 10$ and$s_2 < 0.5, s_1$ . Otherwise the interval ispre_founder(linear over the tested range, slope$t > 3$ ),post_founder(no significant increase of ROP with WOB) orinsufficient_wob_range(5th to 95th percentile WOB span below 4 klbf). - Gives a 90 % residual-bootstrap interval for
$W_b$ . Its resolution is limited by the breakpoint grid. - Recommends
$0.9 \times W_b$ . If the interval is still linear, it recommends the highest tested WOB. It never extrapolates.
The F threshold, slope ratio and 0.9 factor are this project's choices, not values taken from the literature.
RoadMapOptimizer.optimize_by_formation runs one analysis per formation. If optimize_interval gets an interval with several formations it raises an error. For each formation:
- WOB setpoint: the founder recommendation, clipped to the tested range and rig limits, and reduced if a linear torque-vs-WOB fit projects more than 95 % of the top-drive rating.
-
Projected ROP: $\text{ROP}{\text{base}} + [\hat{f}(\text{WOB}{\text{rec}}) - \hat{f}(\text{WOB}_{\text{base}})]$, where
$\hat{f}$ is the fitted ROP-WOB model. It can be below baseline. With no WOB exploration there is no projection and the gain is 0. - Time and cost: hours = length / ROP and cost = hours × day rate / 24. Both can be negative.
- RPM and flow are not optimized. There is no ROP-RPM model and no hydraulics model behind them. If a measured SSI from a time window is passed in and exceeds the limit, the result carries an advisory, but no projected SSI is computed.
These projections are in-sample, made from the same data used for the fit.
For a uniform string fixed at surface and free at the bit, the quarter-wave torsional frequency is
This ignores the BHA. A heavy BHA on a long string behaves closer to a torsional pendulum,
Stick-Slip Severity Index, computed over a time window of downhole (bit) RPM:
This repository classifies SSI below 0.5 as Mild, 0.5 to 1.0 as Moderate and 1.0 or more as Severe. These bands are this project's own thresholds, not an IADC standard. SSI = 1 means the bit momentarily stops when the oscillation is symmetric about the mean. With a speed-controlled top drive, surface RPM stays almost constant during stick-slip. The oscillation then appears in surface torque, so SSI from surface RPM is not meaningful.
Stick-slip periods are a few seconds. A depth log sampled every 0.2 m holds one sample every 30 to 60 s while drilling, so it cannot resolve stick-slip. analyze_stick_slip_window takes a time column in seconds, infers the sampling rate from it, and rejects non-uniform axes. The generator produces such windows with a one-degree-of-freedom torsional pendulum: drillpipe spring, lumped BHA inertia, and bit friction that weakens with speed.
(Jorden and Shirley, 1966; Rehm and McClendon, 1971.) Departures below the normal compaction trend of
-
$\text{TFA} = \sum \frac{\pi}{4}\left(\frac{d_i}{32}\right)^2$ [in²] -
$v_n = \frac{0.32086 \cdot Q}{\text{TFA}}$ [ft/s] -
$\Delta P_{\text{bit}} = \frac{\rho \cdot Q^2}{12032 \cdot C_d^2 \cdot \text{TFA}^2}$ [psi] -
$\text{HSI} = \frac{\Delta P_{\text{bit}} \cdot Q}{1714 \cdot A_b}$ [hp/in²] -
$F_j = \frac{\rho \cdot Q \cdot v_n}{1930.4}$ [lbf]
The generator encodes the founder as the end of the linear ROP-WOB segment. Above it, the incremental slope drops to 20 % and the response rolls over. The analyzer looks for exactly that shape, so the table below shows that the estimator works on data that follow its assumptions, with noise, a sawtooth WOB sweep and finite sampling. It does not validate the physics.
The "quadratic peak" column is the WOB at the maximum of a parabola fitted to the same data, a common simpler definition, shown for comparison.
Seed 42, the committed benchmark (python -c "from src.optimizer.validation import founder_recovery_table; print(founder_recovery_table(42))"):
| Formation | True founder (klbf) | Estimated (klbf) | 90 % interval | Error | Quadratic peak error | Projected ROP change | Noise-free ROP change |
|---|---|---|---|---|---|---|---|
| Shale | 24.0 | 24.25 | 24.25 – 24.25 | +1.0 % | +23.1 % | +7.0 % | +6.8 % |
| Sandstone | 32.0 | 32.21 | 32.05 – 32.37 | +0.7 % | +33.8 % | +6.4 % | +6.6 % |
| Limestone | 38.0 | 38.13 | 37.68 – 38.58 | +0.3 % | +47.5 % | +18.8 % | +18.5 % |
| Dolomite | 42.0 | 42.18 | 41.49 – 42.88 | +0.4 % | +61.6 % | +14.9 % | +15.1 % |
| Sandstone_Lower | 30.0 | 30.30 | 30.14 – 30.47 | +1.0 % | +33.2 % | +18.0 % | +18.0 % |
Over 20 noise seeds (founder_recovery_study(range(20))), the founder was detected in all 100 formation runs. The mean absolute error was 0.8 % to 1.2 % per formation and the worst case was 3.0 % (dolomite). The estimates are biased high by about 1 % on average. The shale interval collapses to a single grid point because the noise there is small relative to the signal. The ROP gains above are measured against the median of an artificial WOB sweep, not against an operating baseline, so they are not expected savings.
Earlier versions of this repository showed a limestone/dolomite misclassification. It came from simulated stick-slip being injected into the surface RPM of the depth log, which distorted the ROP-WOB response. The depth log now has steady surface RPM, and stick-slip is simulated separately in time windows.
Dolomite parameters (WOB 38 klbf, friction 0.52), 2,000 m string, 60 s at 50 Hz, seed 42:
| Surface RPM | SSI (bit RPM) | Class | Dominant period | Time stalled |
|---|---|---|---|---|
| 80 | 1.43 | Severe | 4.0 s | 32 % |
| 120 | 0.20 | Mild | 4.0 s | 0 % |
| 160 | 0.01 | Mild | 4.0 s | 0 % |
The quarter-wave estimate for 2,000 m is
drilling-optimization-engine/
├── data/
│ ├── raw/ # Field data goes here (git-ignored)
│ └── processed/ # Synthetic benchmark: Parquet and LAS 2.0 (seed 42)
├── src/
│ ├── data/
│ │ ├── unit_converter.py # Oilfield / SI unit conversions
│ │ ├── loader.py # LAS (with unit conversion), Parquet and CSV ingestion
│ │ ├── rig_state.py # Rule-based rig state classification
│ │ └── generator.py # Synthetic depth log and torsional time windows
│ ├── physics/
│ │ ├── mse.py # Teale MSE, Dupriest adjustment, efficiency
│ │ ├── d_exponent.py # d, dc-exponent, compaction trend, Eaton
│ │ └── hydraulics.py # TFA, nozzle velocity, bit ΔP, HSI, jet impact
│ ├── dynamics/
│ │ ├── stick_slip.py # SSI, FFT / Welch, f0, time-window analysis
│ │ ├── bit_bounce.py # WOB variability proxies (see Limitations)
│ │ └── wear_model.py # Relative excess-energy wear index
│ └── optimizer/
│ ├── founder_point.py # Hinge-regression founder detection
│ ├── road_map.py # Per-formation WOB setpoint and ROP projection
│ └── validation.py # Founder recovery check on the benchmark
└── tests/ # pytest suite
The loader maps common mnemonics to canonical channels. For LAS files it reads each curve's unit from the header and converts to the canonical unit. A blank or unrecognized unit raises UnitError; it can be declared with unit_overrides={"wob": "klbf"}. CSV and Parquet carry no unit metadata, so they are assumed canonical unless units={...} is passed. TVD is never used as depth.
| Canonical channel | Canonical unit | Also accepted (converted) | Recognized mnemonics |
|---|---|---|---|
depth |
m | ft | DEPTH, DEPT, MD, MEASURED_DEPTH, BIT_DEPTH, HOLE_DEPTH |
rop |
m/hr | ft/hr | ROP, ROPA, ROP_AVG, RATE_OF_PENETRATION, ROP5, ROPM |
wob |
klbf | lbf, kN, N, tonne | WOB, SWOB, DWOB, WEIGHT_ON_BIT, WOB_AVG, WOB_SURF |
rpm |
rev/min | — | RPM, SRPM, DRPM, SURF_RPM, BIT_RPM, ROTARY_SPEED |
torque |
ft-lbf | kft-lbf, N·m, kN·m | TORQUE, TORQ, STOR, SURF_TORQUE, BIT_TORQUE, TORQUE_AVG |
flow |
gpm | L/min, m³/min, m³/s | FLOW, FLOW_IN, FLOWIN, PUMP_FLOW, FLOW_RATE, Q_PUMP |
spp |
psi | bar, kPa, MPa | SPP, STANDPIPE_PRESSURE, PUMP_PRESS, PUMP_PRESSURE, PRESS_SPP |
gamma |
gAPI | — | GR, GAM, GAMMA, GR_CAL, GAMMA_RAY, GRD |
The exact accepted spellings are listed in UNIT_CONVERSIONS in src/data/loader.py.
Requires Python 3.11 or higher.
git clone https://github.com/meilantiago/drilling-optimization-engine.git
cd drilling-optimization-engine
python -m venv venv
source venv/bin/activate # Windows: .\venv\Scripts\Activate.ps1
pip install -r requirements.txtfrom src.data.loader import DrillingDataLoader
from src.optimizer.road_map import RoadMapOptimizer
df = DrillingDataLoader.load_file("data/processed/benchmark_drilling.parquet")
for formation, plan in RoadMapOptimizer().optimize_by_formation(df).items():
print(plan.summary)Regenerate the benchmark files (same seed, identical output):
python -c "from src.data.generator import generate_benchmark_files; generate_benchmark_files()"Run the founder recovery check:
python -m src.optimizer.validationpython -m pytest -qThe suite has 70 tests. On the author's machine it runs in about 3 seconds. The suite covers:
- Unit and formula checks.
- Pandas Series inputs.
- LAS unit conversion and the errors raised for blank or unknown units.
- Rig-state consistency between the per-row and per-DataFrame paths.
- Founder recovery on a known curve and per lithology.
- ROP projections that can be negative, and zero gain without WOB exploration.
- Stick-slip on time windows.
- An end-to-end run.
- Synthetic data only. The benchmark well and the torsional windows come from this repository's generator. No field data has been run through the pipeline.
- Founder validation is partly circular. The generator builds the founder with the same shape the analyzer looks for. The recovery errors above measure the estimator under its own assumptions, not its accuracy on real rock. Real ROP-WOB data are noisier and also respond to bit wear, RPM, hydraulics and lithology changes within an interval.
- ROP projections are in-sample. They come from the fitted ROP-WOB response of the same interval. They cover only WOB changes within the tested range, at the interval's median RPM and flow. On the benchmark, the baseline is the median of an artificial WOB sweep.
- RPM, flow and hydraulics are not optimized. There is no ROP-RPM model, no parasitic pressure losses, no ECD and no rheology. SPP is not used as a constraint.
- Stick-slip needs time data. SSI and the spectra are meaningful only on uniformly sampled time windows of downhole RPM or surface torque. The pendulum model in the generator is illustrative, not calibrated to any rig.
- Bit bounce functions are proxies. On surface WOB from a depth log,
bit_bounce.pymeasures load variability between samples, not axial vibration. Axial vibration needs downhole high-frequency data. - Bit wear is a relative index. It is excess energy over UCS scaled to 0-8 with an uncalibrated capacity. It is not an IADC dull grade. MSE from surface torque includes drag, and UCS is used where confined strength would be more appropriate.
- Rig state uses WOB, RPM, flow and ROP thresholds only. It does not compare bit depth with hole depth.
- d-exponent uses static mud weight (not ECD) and a linear normal-compaction trend, and the generator does not model overpressure.
- Teale, R. (1965). The Concept of Specific Energy in Rock Drilling. International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts, 2(1), 57-73.
- Dupriest, F. E., & Koederitz, W. L. (2005). Maximizing Drill Rates with Real-Time Surveillance of Mechanical Specific Energy. SPE/IADC Drilling Conference, SPE/IADC 92194.
- Jorden, J. R., & Shirley, O. J. (1966). Application of Drilling Performance Data to Overburden and Formation Pressure Detection. Journal of Petroleum Technology, 18(11). SPE 1407.
- Rehm, W. A., & McClendon, R. (1971). Measurement of Formation Pressure From Drilling Data. SPE Annual Fall Meeting, SPE 3601.
- Eaton, B. A. (1975). The Equation for Geopressure Prediction from Well Logs. SPE Annual Fall Meeting, SPE 5544.
- Brett, J. F. (1992). The Genesis of Bit-Induced Torsional Drillstring Vibrations. SPE Drilling Engineering, 7(3), 168-174.
Tiago Meilan
- Degree: Licenciatura en Ciencias Físicas, Universidad de Buenos Aires (UBA).
- Fields of interest: Applied Probability & Statistics, Mathematical Physics, Drilling Mechanics & Geomechanics.
- GitHub: github.com/meilantiago
- LinkedIn: linkedin.com/in/tiago-meilan
- Email: meilantiago@gmail.com