# Math specification — literature backbone

**Status:** Living document. Every implemented term should trace to a citation below or to an explicit `heuristic` calibration note in code + `docs/engineering/validation.md`.

**Primary references:** see `references/papers/README.md` for citation details and source links.

- Romero, A., and Almeida, A. R. (2014). *A Numerical Sucker-Rod Pumping Analysis Tool*. SPE Artificial Lift Conference - Latin America and the Caribbean. [SPE-169395-MS](https://doi.org/10.2118/169395-MS)
- Everitt, T. A., and Jennings, J. W. (1992). *An Improved Finite-Difference Calculation of Downhole Dynamometer Cards for Sucker-Rod Pumps*. SPE Production Engineering. [SPE-18189-PA](https://doi.org/10.2118/18189-PA)
- Araujo, O., et al. (SPE-173970). *3D Rod String Dynamics in Deviated Wells* (minimum-curvature and coupled dynamics reference used for trajectory coupling).
- Eisner, B., Langbauer, C., and Fruhwirth, R. (2022). *A finite element approach for dynamic sucker-rod diagnostics* (Newmark + Rayleigh + diagnostic iteration basis).
- Lukasiewicz, H. (as summarized in Eisner et al.) coupled axial/lateral force balance for deviated rod strings.

---

## 1. Gibbs one-dimensional damped wave (vertical / uniform rod)

**Romero & Almeida — Eq. (2)** (viscous damped wave, constant \(A,E,\rho\)):

\[
\frac{\partial^2 u}{\partial t^2} = a^2 \frac{\partial^2 u}{\partial x^2} - \varsigma \frac{\partial u}{\partial t}
\]

with \(a^2 = E/\rho\) (wave speed) and \(\varsigma\) viscous damping coefficient per unit mass (paper uses lumped fluid damping narrative).

**Variable cross-section (Everitt & Jennings — Eq. 2 form)** after multiplying through by \(\rho A\); axial stiffness gradient:

\[
\frac{\partial}{\partial x}\left( EA \frac{\partial u}{\partial x} \right) = \rho A \frac{\partial^2 u}{\partial t^2} + c \rho A \frac{\partial u}{\partial t} - f_{\text{body}}
\]

where \(f_{\text{body}}\) collects distributed **weight and buoyancy** along the rod (gravity along tangent; buoyancy from fluid — **Lukasiewicz** axial force balance in Eisner Eq. (2) discussion / Lukasiewicz coupled model referenced in Eisner).

**Implementation note:** Code uses discrete \(E_i, A_i, \rho_i\) per node or segment, harmonic or measured surface \(u(0,t)\), and bottom boundary coupling (pump / valve / spring-damper per roadmap).

---

## 2. Damping conventions

- **Gibbs dimensionless damping** \(\nu\): related to decay rate; Romero cites \(\varsigma\) proportional to velocity; Eisner Eq. (1) uses \(\frac{\pi a \nu}{2L}\frac{\partial u}{\partial t}\) form for vertical damped wave.
- **Rayleigh damping (FEA, Eisner):** \(\mathbf{C} = \alpha \mathbf{M} + \beta \mathbf{K}\) for bar elements.
- **XML factors:** `UpStrokeDampingFactor`, `DownStrokeDampingFactor`, `NonDimensionalFluidDamping` modulate effective \(\gamma\) or \(\alpha,\beta\) per phase (mapped in solver; see code comments).

---

## 3. Explicit finite-difference stencil (diagnostic / deviated extension)

**Everitt & Jennings — Eq. (3)** (conceptual explicit FD recursion transferring displacements downhole from surface card). The paper gives a five-point relation in \((i,j)\) space (space index \(i\), time \(j\)) with coefficients involving \(a\), \(c\), \(\Delta t\), \(\Delta x\), and \(EA\).

**Variable \(EA\):** harmonic mean or segment fluxes:

\[
\text{flux}_{i+\frac12} = E_{i+\frac12} A_{i+\frac12} \left( u_{i+1} - u_i \right)
\]

Laplacian-like term formed from \(\partial/\partial x (EA \partial u/\partial x)\) discretization (matches current diagnostic JS prior to C port).

**CFL:** explicit wave requires \(\Delta t \le \text{CFL} \cdot \Delta x / a_{\max}\) with \(a_{\max} = \max \sqrt{E/\rho}\). Code clamps effective CFL (see `verbose.numerics.cflEffective`).

---

## 4. Deviated wells — coupled axial / lateral (reference)

**Lukasiewicz (via Eisner Eqs. (2)–(3))** — force balance along rod tangent and normal; includes \(\rho g A \cos\phi\), viscosity term \(\nu\), Coulomb \(\mu N\), curvature \(R\), and lateral equilibrium.

**Current implementation:** distributed **side load** \(N_i\) and inclination-aware **Coulomb friction** \(F_{f,i} = \mu_{\text{eff}} N_i \operatorname{sgn}(v_i)\) are computed per node and injected into FDM/diagnostic/FEA updates. API exposes `profiles.sideLoadProfile` and `profiles.frictionProfile` when `options.enableProfiles=true`.

---

## 5. Three-dimensional trajectory (SPE-173970)

**Araujo et al. — Eqs. (5)–(24):** minimum-curvature method between survey stations: unit tangent, curvature angle \(\gamma_i\), radius \(r_{c,i}\), binormal, center of curvature, position propagation \(R(s_{i+1})\), and interpolation of any MD \(s\) within a segment.

**Implementation:** `solver-c/src/trajectory.c` uses minimum-curvature style tangent interpolation and exports node-wise \(\kappa(s)\), inclination, and azimuth on the rod grid for side-load/friction coupling and `trajectory3D` output.

---

## 6. Dynamic 1D bar FEM (Eisner et al.)

- **Bar stiffness / mass:** consistent element \(K_e\), \(M_e\) for linear shape functions.
- **Newmark-β** (\(\beta=\frac14\), \(\gamma=\frac12\)) for transient axial dynamics.
- **Bottom BC:** spring-damper + friction; diagnostic mode adds **iterative plunger load** adjustment so that computed top reaction matches measured polished load (Eisner Fig. 4 principle — restart / bisection per time step).

---

## 7. Boundary conditions

| Mode | Surface | Bottom |
|------|---------|--------|
| **Predictive** | Harmonic \(u(0,t)\) (Romero **Eq. (4)** style crank motion) or measured surrogate | Pump: spring-damper + friction + valve-state pressure balance; gas fraction inferred from chamber state |
| **Diagnostic** | Measured \(u(0,t)\) and \(F(0,t)\) from dynamometer card (Everitt) | Same pump BC; FEA adjusts unknown plunger load and reports valve/gas timeline |

**Surface card QA (Eisner narrative):** ≥ 75 samples recommended for Fourier-style tools; we enforce minimum samples and cycle checks in API (`POST /solve/validate-card`).

---

## 8. Outputs (API / solver)

- Polished and downhole **cards** \((x, F)\).
- **Pump movement:** plunger position series, velocity, stroke, period.
- **Optional profiles:** `stressProfile`, `sideLoadProfile`, `trajectory3D` (when `schemaVersion >= 2`).
- **Comparison:** FDM vs FEA peak deltas + point-wise residuals + optional Fourier analytical baseline (`comparison.fourier` when enabled).

---

## 9. Symbols (SI internal)

| Symbol | Unit | Meaning |
|--------|------|---------|
| \(u\) | m | Axial displacement |
| \(F\) | N | Axial force (tension positive) |
| \(E\) | Pa | Young’s modulus |
| \(A\) | m² | Cross-sectional area |
| \(\rho\) | kg/m³ | Rod density |
| \(\phi\) | rad | Inclination from vertical |
| \(\kappa\) | 1/m | Path curvature |

---

## 10. Implementation quality rules (for next execution)

- Every newly introduced solver term must be tagged as one of:
  - `equation-backed` (with paper/equation reference),
  - `heuristic` (with rationale + planned replacement).
- No silent heuristic drift: if a coefficient changes, update validation fixtures and notes.
- Any new field wired into equations must be reflected in `docs/engineering/field-traceability.md`.
- Any new coupled term must include at least:
  - one unit-scale numerical test,
  - one integration case assertion,
  - one regression guard.

---

*End of math spec backbone. Extend with equation numbers from PDFs as features land.*