Bit Run Analysis - Footage, ROP, MSE & Cost per Foot

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 pandas as pd


# ── Verified Mechanical Specific Energy (do not edit) ────────────────────
def compute_mse(torque_ftlbs, rpm, wob_lbf, rop_fthr, bit_diameter_in):
    """
    Calculate Mechanical Specific Energy.

    Parameters
    ----------
    torque_ftlbs : array-like
        Torque (ft-lbs).
    rpm : array-like
        Rotary speed (rev/min).
    wob_lbf : array-like
        Weight on bit (lbf). Note: input in lbf, not klbf.
    rop_fthr : array-like
        Rate of penetration (ft/hr).
    bit_diameter_in : float
        Bit diameter (inches).

    Returns
    -------
    numpy.ndarray
        MSE values (psi).
    """
    torque = np.asarray(torque_ftlbs, dtype=float)
    r = np.asarray(rpm, dtype=float)
    w = np.asarray(wob_lbf, dtype=float)
    rop = np.asarray(rop_fthr, dtype=float)
    db = bit_diameter_in

    # Avoid division by zero
    rop_safe = np.where(rop > 0, rop, 1e-6)

    rotary_term = 480 * torque * r / (db**2 * rop_safe)
    wob_term = 4 * w / (np.pi * db**2)

    return rotary_term + wob_term


# ── Drilling economics constants (do not edit) ───────────────────────────
BIT_DIAMETER = 8.5               # bit diameter, in
TIME_INTERVAL_HR = 2 / 60        # 2-minute log sample interval, hr
SPREAD_RATE_PER_HR = 280000 / 24  # $280k/day rig spread rate, USD/hr

BIT_RUNS = [(0, 3500), (3500, 7200), (7200, 9800)]  # depth windows, ft
BIT_COSTS = [22000.0, 28000.0, 35000.0]             # bit invoice cost, USD


# ── Self-contained OML-58 well log (do not edit) ─────────────────────────
def _build_well():
    """Build a deterministic ~9,800-ft drilling log (chapter-style data gen).

    The log is sampled at a uniform 2-minute interval and Depth is the TRUE
    cumulative sum of ROP * TIME_INTERVAL_HR (no rescaling), so the drilling
    identity footage = ROP * hours holds exactly and the well reaches
    ~9,800 ft on its own. The formation hardens with depth: a soft fast top,
    a medium middle, and a hard slow base -- so ROP falls and torque/MSE
    climb downhole.
    """
    np.random.seed(143)
    n = 5350
    x = np.linspace(0, 1, n)
    rop_base = np.piecewise(
        x,
        [x < 0.3, (x >= 0.3) & (x < 0.6), x >= 0.6],
        [lambda t: 120 - 80 * t,           # soft formation, fast drilling
         lambda t: 45 + 10 * np.sin(20 * t),  # medium formation, variable
         lambda t: 25 - 10 * (t - 0.6)],   # hard formation, slow
    )
    rop = rop_base + np.random.normal(0, 5, n)
    rop = np.clip(rop, 5, 200)

    # True cumulative depth from ROP at the 2-minute interval (no normalization)
    # so footage = ROP * hours stays physically consistent across each bit run.
    depth = np.cumsum(rop * TIME_INTERVAL_HR)

    wob = np.where(depth < 4000, 15 + np.random.normal(0, 2, n),
          np.where(depth < 7000, 25 + np.random.normal(0, 3, n),
                   35 + np.random.normal(0, 3, n)))
    wob = np.clip(wob, 5, 50)

    rpm = np.where(depth < 4000, 140 + np.random.normal(0, 10, n),
          np.where(depth < 7000, 120 + np.random.normal(0, 8, n),
                   100 + np.random.normal(0, 8, n)))
    rpm = np.clip(rpm, 40, 180)

    torque = 0.125 * (wob * 1000.0) * (BIT_DIAMETER / 12.0) + np.random.normal(0, 250, n)
    torque = np.clip(torque, 800, 30000)

    return pd.DataFrame({
        "Depth": np.round(depth, 1),
        "ROP": np.round(rop, 1),
        "WOB": np.round(wob, 1),
        "RPM": np.round(rpm, 0).astype(int),
        "Torque": np.round(torque, 0).astype(int),
    })


well = _build_well()


def analyze_bit_run(df, depth_lo, depth_hi, bit_cost,
                    db=BIT_DIAMETER, interval_hr=TIME_INTERVAL_HR,
                    spread_rate=SPREAD_RATE_PER_HR):
    """Summarize one bit run over the depth window [depth_lo, depth_hi)."""
    sub = df[(df.Depth >= depth_lo) & (df.Depth < depth_hi)]
    footage = sub.Depth.max() - sub.Depth.min()
    hours = len(sub) * interval_hr
    avg_rop = sub.ROP.mean()
    avg_mse = compute_mse(sub.Torque, sub.RPM, sub.WOB * 1000, sub.ROP, db).mean()
    cost_per_foot = (bit_cost + spread_rate * hours) / footage
    return {
        "footage": float(footage),
        "hours": float(hours),
        "avg_rop": float(avg_rop),
        "avg_mse": float(avg_mse),
        "cost_per_foot": float(cost_per_foot),
    }


def best_bit(df, runs, costs, db=BIT_DIAMETER, interval_hr=TIME_INTERVAL_HR,
             spread_rate=SPREAD_RATE_PER_HR):
    """Return the integer index of the run with the lowest cost per foot."""
    cpfs = [analyze_bit_run(df, lo, hi, cost, db, interval_hr, spread_rate)["cost_per_foot"]
            for (lo, hi), cost in zip(runs, costs)]
    return int(np.argmin(cpfs))


run_results = [analyze_bit_run(well, lo, hi, cost)
               for (lo, hi), cost in zip(BIT_RUNS, BIT_COSTS)]
best_idx = best_bit(well, BIT_RUNS, BIT_COSTS)

for i, r in enumerate(run_results):
    print(f"Bit {i + 1}: footage={r['footage']:.1f} ft  hours={r['hours']:.2f}  "
          f"avg_ROP={r['avg_rop']:.1f} ft/hr  avg_MSE={r['avg_mse']:,.0f} psi  "
          f"$/ft={r['cost_per_foot']:.2f}")
print("Best bit (lowest $/ft):", best_idx)

lockCopying code is a Full Access feature.