Exerciseschevron_rightChapter 15chevron_right15.2
fitness_center

Exercise 15.2

Gas Property Module - Z, Bg, Density, Viscosity & Compressibility Table

Level 2
Chapter 15: Gas Engineering
descriptionProblem

Package the chapter's property functions into a single gas_properties(P, T, gamma_g) function that returns a dict of {Z, Bg, density, viscosity, cg} in one call, then add a gas_property_table(P_range, T, gamma_g) wrapper that returns a formatted DataFrame. Validate against published data for a 0.7 gravity gas at 250°F.

---

You'll assemble the chapter's gas-PVT toolkit into a single property table. The six verified correlations are embedded for you; do not modify them or re-derive the physics. The Sutton pseudo-criticals, the Dranchuk–Abou-Kassem Z-factor, gas_fvf, gas_density, the Lee–Gonzalez–Eakin viscosity, and the numerical gas_compressibility are all in scope already.

Your job is the one missing piece: the table builder.

Write gas_property_table(P_range, T_F, gamma_g) that returns a pandas.DataFrame with one row per pressure in P_range and these columns, in this exact spelling and order:

'P (psia)', 'Z', 'Bg (rcf/scf)', 'rho (lb/ft3)', 'mu (cp)', 'cg (1/psi)'

For each pressure P:

  1. Get Ppc, Tpc = pseudo_critical_properties(gamma_g).
  2. Tpr = (T_F + 459.67) / Tpc.
  3. Z = z_factor_DAK(P / Ppc, Tpr).
  4. Bg = gas_fvf(P, T_F, Z) → column 'Bg (rcf/scf)'
  5. rho = gas_density(P, T_F, Z, gamma_g) → column 'rho (lb/ft3)'
  6. mu = gas_viscosity_lee(P, T_F, Z, gamma_g) → column 'mu (cp)'
  7. cg = gas_compressibility(P, T_F, Z, gamma_g) → column 'cg (1/psi)'

Store the raw P in 'P (psia)' and the same Z you computed in 'Z'. Keep the numbers as floats (no rounding/formatting; the grader reads the raw values).

The embedded constants are the chapter's published validation case (a 0.70-gravity gas at 250 °F):

GAMMA_G = 0.70
T_F     = 250.0
P_RANGE = [2000.0, 4000.0, 6000.0]   # psia

Then build the table and expose these exact names the tests read:

props    = gas_property_table(P_RANGE, T_F, GAMMA_G)   # the DataFrame
z_2000   = float(props.loc[0, 'Z'])
bg_2000  = float(props.loc[0, 'Bg (rcf/scf)'])
rho_4000 = float(props.loc[1, 'rho (lb/ft3)'])

> Think about it: as you walk up the pressure range, Bg should fall > (gas compresses into less reservoir volume per scf) while rho rises. The > ideal-gas skeleton Bg = 0.02827 · Z · T_R / P means Bg · P / Z is the same > constant on every row, a quick sanity check that your Z is threaded through > correctly. Why does cg shrink toward roughly 1/P as pressure climbs?

lightbulbHints (0/3)

Stuck? Reveal hints one at a time — they progress from nudge to near-solution.

codeYour solution
main.py
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 pandas as pd
from scipy.optimize import brentq


# ── Verified gas-PVT correlations from Chapter 15 (do not edit) ──────────
def pseudo_critical_properties(gas_sg):
    """
    Estimate pseudo-critical properties from gas specific gravity.
    Uses the Sutton (1985) correlations.

    Parameters
    ----------
    gas_sg : float
        Gas specific gravity (air = 1.0). Typical range: 0.55 to 0.90.

    Returns
    -------
    tuple
        (Ppc in psia, Tpc in deg R)
    """
    Ppc = 756.8 - 131.0 * gas_sg - 3.6 * gas_sg**2
    Tpc = 169.2 + 349.5 * gas_sg - 74.0 * gas_sg**2
    return Ppc, Tpc


def z_factor_DAK(Ppr, Tpr):
    """
    Calculate gas Z-factor using the Dranchuk-Abou-Kassem (1975)
    equation of state.

    This correlation fits the Standing-Katz chart with 11 empirical
    coefficients. Valid for 1.0 <= Tpr <= 3.0 and 0.2 <= Ppr <= 15.

    Parameters
    ----------
    Ppr : float
        Pseudo-reduced pressure (dimensionless).
    Tpr : float
        Pseudo-reduced temperature (dimensionless).

    Returns
    -------
    float
        Compressibility factor Z.
    """
    # DAK (1975) coefficients
    A1, A2, A3 = 0.3265, -1.0700, -0.5339
    A4, A5, A6 = 0.01569, -0.05165, 0.5475
    A7, A8 = -0.7361, 0.1844
    A9, A10, A11 = 0.1056, 0.6134, 0.7210

    def residual(rho_r):
        """DAK residual function. Z = 0.27 * Ppr / (rho_r * Tpr)"""
        Z = 0.27 * Ppr / (rho_r * Tpr)

        F = (
            A1 + A2/Tpr + A3/Tpr**3 + A4/Tpr**4 + A5/Tpr**5
        ) * rho_r + (
            A6 + A7/Tpr + A8/Tpr**2
        ) * rho_r**2 - A9 * (
            A7/Tpr + A8/Tpr**2
        ) * rho_r**5 + A10 * (
            1 + A11*rho_r**2
        ) * (rho_r**2 / Tpr**3) * np.exp(-A11*rho_r**2) + 1 - Z

        return F

    # Solve for reduced density
    rho_r = brentq(residual, 0.001, 5.0)
    Z = 0.27 * Ppr / (rho_r * Tpr)
    return Z


def gas_fvf(P_psia, T_F, Z):
    """Gas formation volume factor (rcf/scf -- reservoir ft^3 per scf).

    The 0.02827 constant gives reservoir CUBIC FEET per scf. For reservoir
    bbl/scf, divide by 5.615 (i.e. use 0.005035).
    """
    T_R = T_F + 459.67
    return 0.02827 * Z * T_R / P_psia


def gas_density(P_psia, T_F, Z, gas_sg):
    """Gas density (lb/ft3)."""
    T_R = T_F + 459.67
    return 2.7 * gas_sg * P_psia / (Z * T_R)


def gas_viscosity_lee(P_psia, T_F, Z, gas_sg):
    """
    Gas viscosity using Lee-Gonzalez-Eakin (1966) correlation.

    Returns viscosity in centipoise (cp).
    """
    T_R = T_F + 459.67
    Mg = 28.97 * gas_sg  # molecular weight
    rho_g = gas_density(P_psia, T_F, Z, gas_sg)  # lb/ft3
    rho_gcc = rho_g / 62.428  # convert to g/cc

    K = (9.4 + 0.02 * Mg) * T_R**1.5 / (209 + 19*Mg + T_R)
    X = 3.5 + 986/T_R + 0.01*Mg
    Y = 2.4 - 0.2*X

    return 1e-4 * K * np.exp(X * rho_gcc**Y)


def gas_compressibility(P_psia, T_F, Z, gas_sg, dP=10):
    """
    Gas isothermal compressibility cg (1/psi).
    Calculated numerically from Z(P).
    """
    Ppc, Tpc = pseudo_critical_properties(gas_sg)
    T_R = T_F + 459.67
    Tpr = T_R / Tpc

    Ppr_plus = (P_psia + dP) / Ppc
    Ppr_minus = (P_psia - dP) / Ppc

    try:
        Z_plus = z_factor_DAK(Ppr_plus, Tpr)
        Z_minus = z_factor_DAK(Ppr_minus, Tpr)
        dZ_dP = (Z_plus - Z_minus) / (2 * dP)
        cg = 1/P_psia - 1/Z * dZ_dP
    except (ValueError, RuntimeError):
        cg = 1/P_psia  # ideal gas approximation as fallback

    return cg


# ── Validation constants: 0.70-gravity gas at 250 deg F (do not edit) ────
GAMMA_G = 0.70             # gas specific gravity (air = 1.0)
T_F = 250.0               # reservoir temperature, deg F
P_RANGE = [2000.0, 4000.0, 6000.0]   # pressures, psia


def gas_property_table(P_range, T_F, gamma_g):
    """Build a gas-property table, one row per pressure.

    Columns (exact spelling/order):
        'P (psia)', 'Z', 'Bg (rcf/scf)', 'rho (lb/ft3)', 'mu (cp)', 'cg (1/psi)'
    """
    Ppc, Tpc = pseudo_critical_properties(gamma_g)
    T_R = T_F + 459.67
    Tpr = T_R / Tpc

    rows = []
    for P in P_range:
        Z = z_factor_DAK(P / Ppc, Tpr)
        rows.append({
            "P (psia)": P,
            "Z": Z,
            "Bg (rcf/scf)": gas_fvf(P, T_F, Z),
            "rho (lb/ft3)": gas_density(P, T_F, Z, gamma_g),
            "mu (cp)": gas_viscosity_lee(P, T_F, Z, gamma_g),
            "cg (1/psi)": gas_compressibility(P, T_F, Z, gamma_g),
        })
    return pd.DataFrame(rows)


props = gas_property_table(P_RANGE, T_F, GAMMA_G)
z_2000 = float(props.loc[0, "Z"])
bg_2000 = float(props.loc[0, "Bg (rcf/scf)"])
rho_4000 = float(props.loc[1, "rho (lb/ft3)"])

print(f"GAS PROPERTY TABLE - gamma_g = {GAMMA_G}, T = {T_F} deg F")
print(props.to_string(index=False))
print()
print("Z @ 2000 psia   :", z_2000)
print("Bg @ 2000 psia  :", bg_2000, "rcf/scf")
print("rho @ 4000 psia :", rho_4000, "lb/ft3")

lockCopying code is a Full Access feature.