Skip to content

MEIO — EOQ Reorder-Point Refinement & the Type-2 Loss Function

Chapter summary

This chapter analyses how the Economic Order Quantity (EOQ) refines the reorder point (ROP) / safety stock inside the MEIO engine, comparing the V1 Python, V2/V3 Rust, and V3 standalone Python implementations. It establishes why the reorder point must decrease as the EOQ increases (a Type-2 / fill-rate effect, not a Type-1 / cycle-service-level effect), documents where each version gets this right and where it is inconsistent or approximate, and presents an improvement plan that replaces the current EOQ averaging heuristic with the exact first-order loss-function model β = 1 − n₁(s)/Q. The companion implementation plan lives at .kilo/plans/1782658808592-meio-eoq-rop-loss-function.md.

Related chapters

  • MEIO Theory — the greedy marginal-value algorithm, BOM/kit wait times, reverse logistics, Erlang-B asset mode.
  • Safety Stock — the fundamental ROP = d̄·L + SS derivation and the Type-1 vs Type-2 distinction.
  • EOQ & Lot-Sizing — how the EOQ itself is computed (Classic EOQ, K-Curve, Wagner-Whitin) and the per-item eoq field.

TOC

  1. The intuition — why ROP must drop as EOQ rises
  2. Type-1 vs Type-2 service — the crux
  3. The (s, Q) loss-function theory
  4. Current state — V1 (Python)
  5. Current state — V2/V3 (Rust)
  6. Current state — V3 standalone Python
  7. Version-comparison matrix
  8. The three problems
  9. The improvement plan
  10. Risks & validation
  11. Key files

The intuition — why ROP must drop as EOQ rises

Consider a continuous-review (s, Q) policy: when the inventory position drops to the reorder point s, an order of size Q (the EOQ) is placed. The order arrives after the lead time L. During that lead time, demand continues to arrive and the on-hand inventory is consumed.

The key insight is that a larger order quantity spreads the risk of a stockout over a longer cycle. Each replenishment cycle places the system back at inventory position = s + Q. A stockout within the cycle can only occur during the lead-time window, but the fraction of demand that is back-ordered depends on how much demand the cycle must satisfy relative to the shortage:

\[ \text{fill rate} \;\beta = 1 - \frac{\text{expected shortage per cycle}}{\text{demand per cycle}} = 1 - \frac{n_1(s)}{Q} \]

where n₁(s) = E[(X − s)⁺] is the first-order loss function (expected shortage) for demand-during-lead-time X at reorder point s.

Because n₁(s)/Q shrinks as Q grows, β rises with Q at fixed s. Therefore, to hold a target fill rate β* constant, s (the ROP, and hence the safety stock SS = s − d̄L) must decrease as the EOQ increases.

This is the planner's intuition, and it is the textbook result (Silver, Pyke & Peterson; Axsäter). It is also the reason the MEIO engine accepts an eoq field per SKU — the EOQ is not only a cycle-stock input, it is a safety-stock input.


Type-1 vs Type-2 service — the crux

There are two distinct "service level" definitions, and only one of them depends on the EOQ:

Metric Definition EOQ-dependent?
Cycle Service Level (Type-1, CSL) \(P(X \le s) = F(s)\) — probability of no stockout in a cycle No — depends only on s
Fill Rate (Type-2, β) \(1 - n_1(s)/Q\) — fraction of demand satisfied directly from stock Yes — depends on s and Q
flowchart LR
    subgraph Cycle["One replenishment cycle"]
        direction LR
        A["inv pos = s + Q<br/>(order arrives)"] --> B["demand consumes stock<br/>over the cycle"]
        B --> C["inv pos = s<br/>(reorder triggered)"]
        C --> D["lead-time window<br/>demand X ~ D(L)"]
        D --> E{"X ≤ s ?"}
        E -->|yes| F["no stockout<br/>cycle OK"]
        E -->|no| G["shortage = X − s<br/>back-ordered demand"]
    end
    G -.->|"n₁(s) = E[(X−s)⁺]<br/>β = 1 − n₁(s)/Q"| H["Type-2 fill rate"]
    F -.->|"P(X ≤ s) = F(s)"| I["Type-1 cycle SL"]

A common error is to size safety stock using the Type-1 formula SS = z_α · σ · √L (which is EOQ-independent) and then claim the EOQ has no effect on the buffer. That is only true for the cycle service level; the fill rate — which is what MEIO's group targets actually constrain — is genuinely EOQ-sensitive. The MEIO engine's group fill_rate_target is a Type-2 objective, so the EOQ must enter the ROP math.

Worked example

Demand-during-lead-time X ~ Normal(μ=100, σ=30), target fill rate β* = 0.99:

EOQ Q n₁(s) needed = Q(1−β*) s solving n₁(s) = Q(1−β*) Safety stock s − μ
50 0.5 ≈ 100.4 (nearly at the mean) +0.4
100 1.0 ≈ 92.3 −7.7
200 2.0 ≈ 82.6 −17.4
500 5.0 ≈ 70.0 −30.0

As Q grows, the same target fill rate is achieved with a lower reorder point — eventually dropping below the mean lead-time demand (negative safety stock), which is correct: a large enough order quantity back-fills the shortage from the cycle's own surplus.


The (s, Q) loss-function theory

The first-order loss function is the expected shortage above a threshold:

\[ n_1(s) = E[(X - s)^+] = \int_s^\infty (1 - F(x))\, dx \quad \text{(continuous)} \]
\[ n_1(s) = \sum_{k > s} (k - s)\, P(X = k) \quad \text{(discrete)} \]

The exact Type-2 fill rate is then:

\[ \boxed{\;\beta(s, Q) = 1 - \frac{n_1(s)}{Q}, \quad Q > 0\;} \]

and the inverse (find s for a target β*) solves n₁(s) = Q(1 − β*).

Per-distribution loss functions

Distribution of X n₁(s)
Normal N(μ, σ²) σ · L(z), where z = (s−μ)/σ and L(z) = φ(z) − z(1 − Φ(z)) is the standard normal loss
Poisson Pois(λ) discrete sum Σ_{k>s} (k−s)·PMF(k), or recursion n₁(s) = n₁(s−1) − (1 − F(s−1))
Gamma Γ(α, θ), mean μ = αθ μ(1 − F_{α+1}(s)) − s(1 − F_α(s)), where F_k is the Gamma CDF with shape k
LogNormal (μ, σ²) E[X]·Φ((μ+σ²−ln s)/σ) − s·Φ((μ−ln s)/σ), E[X] = exp(μ+σ²/2)
Weibull (k, λ) numerical integral ∫_s^∞ (1 − F(x)) dx (Gauss–Legendre quadrature)
NegativeBinomial (r, p) discrete loss sum / recursion (as Poisson)
Empirical samples {xᵢ} (1/N) Σ max(xᵢ − s, 0) (mean excess)

The Normal case σ·L(z) is the standard "normal loss" tabulated in every inventory textbook; the others are the natural generalisations.


Current state — V1 (Python)

The legacy V1 implementation lives in files/MEIO/intial/. It is not the production MEIO path (production runs the Rust V2/V3 optimizer), but it is the historical reference and is intentionally left unchanged for legacy comparison.

The averaging forward loop

Distributions.py:90golden_formula_ROP2FR and the sibling poisson_distribution_ROP2FR (:32) compute the fill rate by averaging the CDF over the inventory-position range [ROP, ROP + EOQ]:

step = int(max(1, avgsize, math.ceil((eoq)/maxEoqFreq)))   # maxEoqFreq = 10
tested_eoq = 0
while tested_eoq < eoq + step:
    maxEoq = min(eoq, tested_eoq)
    fill_rate = ndtr((safetyStock) / racine)               # normal CDF at ROP+tested_eoq
    total_fill_rate += fill_rate
    tested_eoq += step
return total_fill_rate / loop_counter                       # average over ~10 points

This is a midpoint / trapezoidal approximation of (1/Q) ∫_0^Q F(s + q) dq, which is related to the loss function (the integral form of n₁ is (1/Q) ∫_0^Q (1−F(s+q)) dq) but only coarsely. It effectively subtracts roughly eoq/2 from the buffer.

The −eoq/2 heuristic in the inverse & jump

  • poisson_distribution_FR2ROP (Distributions.py:20): rop = poisson.ppf(FR)·avgsize − int(eoq)/2 — a direct −eoq/2 heuristic.
  • MEIO.initialJump (MEIO.py:131): returned_buffer = max(avg_size, ceil(forecast_ovr_lt·avg_size − ceil(eoq/2))).

The EOQ-blind gaps (V1)

  • normal_distribution_FR2ROP (Distributions.py:78): rop = norm.ppf(FR, loc=forecastlt, scale=lt_stdev)·avgsizeignores EOQ entirely. So high-volume (Normal) SKUs get no EOQ-driven ROP reduction.
  • MEIO.minMaxROPQty (MEIO.py:481): tgt_min_rop_qty = tgt_min_sl_qty + effective_total_lt_fcst — the ROP bounds ignore EOQ, so the EOQ-driven reduction can be capped by the floor.

Current state — V2/V3 (Rust)

The production optimizer lives in files/MEIO_v3/src/ and serves both V2 and V3 (the V3 additions are the empirical/bootstrap distributions; the fill-rate math path is shared). See meio_runner.py:1375 — V3 is the default, V2 is the explicit fallback.

The V2 eoq_step fix

math.rs:175 — the integration step is now driven by the route lot size (effective_lot_size(), i.e. route.mult_qty) instead of the demand-line avg_size. V1 used avg_size, which for weekly demand buckets made the step ~13× too large (meio_runner.py:1164). The forward averaging loop in distributions.rs:300 (golden_formula_rop_to_fr), :363 (poisson_rop_to_fr) and :381 (fitted_cdf_rop_to_fr) is structurally identical to V1 — same averaging approximation, just with a corrected step.

Where EOQ is applied correctly (V2/V3)

  • distributions.rs:402 poisson_fr_to_rop: rop = k·avg_size − eoq/2 — keeps the −eoq/2 heuristic ✓
  • marginal.rs:146 initial_jump: first_buf = max(lot_size, ceil(fcst_lt·lot_size − eoq/2))

The EOQ-blind gaps (V2/V3)

  • distributions.rs:425 normal_fr_to_rop: rop = (forecast_lt + z·σ)·avg_sizeignores EOQ. This is the same defect as V1's normal_distribution_FR2ROP. High-volume Normal SKUs (the majority of the catalogue) get no EOQ-driven ROP reduction.
  • marginal.rs:179 further_jump — the Normal branch (:212) computes the next buffer from the z-score and lead-time-demand only; EOQ is not in the formula. Only the Poisson branch (:201) calls the EOQ-aware poisson_fr_to_rop.
  • marginal.rs:39 min_max_rop_qty: tgt_min_rop_qty = tgt_min_sl_qty + effective_total_lt_fcst — the floor ignores EOQ, so the reduction can be capped.

The V3 empirical path

The V3 empirical / compound-Poisson / bucketed distributions (distributions.rs:42-82) route through fitted_cdf_rop_to_fr (:381), which uses the same averaging loop as the parametric distributions — so the EOQ effect is the same coarse midpoint approximation, not the exact empirical mean-excess loss function.


Current state — V3 standalone Python

files/inventory_engine/rop.py is the V3 reference / dashboard engine (a pure pure-Python empirical-bootstrap inventory optimiser). It is not the production MEIO write path today (the Rust optimizer writes PIPE_meio_results), but it is the documented "future" ROP engine and is exposed via the service-level-curve API (main.py:27560).

rop.py:19reorder_point:

def reorder_point(simulated_ltd_samples, service_level):
    return EmpiricalDistribution(simulated_ltd_samples).inv_cdf(service_level)

rop.py:73fill_rate_for_buffer:

def fill_rate_for_buffer(simulated_ltd_samples, buffer):
    return EmpiricalDistribution(simulated_ltd_samples).cdf(buffer)

Both are pure Type-1 (F(s), percentile / CDF) with no EOQ parameter at all. If this engine is ever wired into the production safety-stock computation without an eoq argument, every SKU's buffer would silently lose the EOQ-driven reduction. The improvement plan (below) adds an eoq=0.0 parameter that preserves Type-1 behaviour when eoq=0 and applies the exact loss function when eoq>0.


Version-comparison matrix

Code path V1 Python V2/V3 Rust V3 Python reference
Forward fill rate fill_rate_for_rop averaging loop (coarse ~−Q/2) averaging loop (V2 lot-size step) CDF(buffer)no EOQ
Inverse ROP (Poisson) ppf(β)·avgsize − Q/2 k·avg_size − Q/2 n/a (percentile)
Inverse ROP (Normal) ignores EOQ ignores EOQ n/a (percentile)
Initial jump max(avg, fcst·avg − Q/2) max(lot, fcst·lot − Q/2) n/a
Further jump (Poisson) ppf(β+Δ) − Q/2 poisson_fr_to_rop n/a
Further jump (Normal) ignores EOQ ignores EOQ n/a
Min/max ROP bounds ignores EOQ ignores EOQ n/a
Empirical / fitted n/a averaging loop (coarse) percentile — no EOQ

✓ = EOQ applied (heuristic). ✗ = EOQ-blind gap.


The three problems

  1. Incompleteness. The Normal inverse (normal_fr_to_rop), the Normal further_jump branch, and the min_max_rop_qty floor are EOQ-blind in every version. Since high-volume SKUs (the bulk of inventory value) use the Normal path, most of the catalogue sees almost no EOQ-driven ROP reduction today.

  2. Approximation quality. Where EOQ is applied, it is via the averaging forward loop (a midpoint proxy that effectively subtracts ~Q/2) and the explicit −eoq/2 heuristic in the inverse. This over-reduces the ROP for large Q relative to the exact loss function 1 − n₁(s)/Q, because the averaging integral (1/Q)∫F(s+q)dq is only the linear part of the loss function and ignores the curvature of n₁ for heavy-tailed distributions.

  3. The V3 percentile engine drops EOQ entirely. inventory_engine/rop.py computes ROP = percentile(demand_over_LT, SL) — a pure Type-1 cycle service level with no Q. If wired into production safety stock, it would silently remove the EOQ benefit.


The improvement plan

Replace the EOQ averaging approximation with the exact Type-2 loss-function fill rate β = 1 − n₁(s)/Q across the Rust MEIO_v3 engine and the V3 Python inventory_engine reference. V1 Python is left as legacy.

Design decisions (resolved with the planner):

  • Scope: Rust files/MEIO_v3/src/ + V3 Python files/inventory_engine/rop.py. V1 Python (files/MEIO/intial/) left untouched for legacy comparison.
  • eoq = 0 (transfer / make routes, no ordering cost): fall back to Type-1 F(s). Preserves current transfer-route buffers (the averaging loop collapses to a single F(s) point when Q=0).
  • Full replacement of the averaging forward loop with the exact loss function — not just an additive fix.

1. New loss module — files/MEIO_v3/src/math.rs (or loss.rs)

Implement n₁(s) and its inverse solve_n1(target, Q, β) → s for every DistributionType (see the per-distribution table above).

2. Forward fill rate — files/MEIO_v3/src/distributions.rs

Rewrite fill_rate_for_rop (:225) to a single exact formula:

if eoq <= 0.0 {
    fr = dist.cdf(buffer)                          // Type-1 fallback (transfer routes)
} else {
    fr = (1.0 - n1(dist, buffer) / eoq).clamp(0.0, 1.0)  // Type-2 exact
}

Delete the tested_eoq averaging loops in golden_formula_rop_to_fr (:300), normal_rop_to_fr (:344), poisson_rop_to_fr (:363), fitted_cdf_rop_to_fr (:381).

3. Inverse ROP — files/MEIO_v3/src/distributions.rs

  • poisson_fr_to_rop (:402): s = solve_n1(…) instead of k·avg_size − Q/2.
  • normal_fr_to_rop (:425): s = solve_n1(…)now EOQ-aware. Key fix for high-volume Normal SKUs.
  • New fitted_fr_to_rop(dist, target, eoq, …) for Gamma/LogNormal/Weibull/NegBin/ Empirical so the fitted & V3-empirical paths also get the exact inverse.

4. Jumps & bounds — files/MEIO_v3/src/marginal.rs

  • initial_jump (:146): derive first_buf from the exact inverse at sku_min_fill_rate (floored at lot_size).
  • further_jump (:179): both branches use the exact inverse — the Normal branch stops ignoring EOQ.
  • min_max_rop_qty (:39): add the EOQ term to the floor so the reduction is not capped: tgt_min_rop_qty = max(0, tgt_min_sl_qty + fcst_lt − eoq/2).

5. V3 Python reference — files/inventory_engine/rop.py

  • reorder_point(samples, service_level, eoq=0.0): when eoq > 0, solve mean_excess(samples, s) = eoq·(1−service_level); else keep the percentile.
  • fill_rate_for_buffer(samples, buffer, eoq=0.0): 1 − mean_excess/eoq when eoq > 0, else empirical CDF.
  • Default eoq=0.0 preserves all existing call sites.

6. Tests & validation

  • Property tests: β monotone in s; β increasing in Q; ROP non-increasing in Q; eoq=0 returns F(s) (Type-1 parity).
  • Numerical parity vs scipy for Normal/Poisson/Gamma loss functions.
  • Dev-pipeline A/B: confirm committed_buffer decreases as EOQ increases and group fill-rate targets still hold.

Full task breakdown: see the implementation plan at .kilo/plans/1782658808592-meio-eoq-rop-loss-function.md.


Risks & validation

Every SKU's buffer changes

The averaging proxy over-reduces the ROP at large Q vs the exact loss function. After the rewrite, most high-EOQ SKUs will see committed_buffer rise toward the exact value (still a decrease vs the eoq=0 baseline — the EOQ benefit remains, just calibrated). Run a dev-pipeline A/B and flag any group that falls below its fill_rate_target.

  • Weibull numerical integral adds cost per fill_rate_for_rop call; the greedy loop calls it many times per SKU. Profile and cache the quadrature nodes.
  • Empirical inverse via mean-excess over sorted samples is O(log N) per evaluation after a one-time sort — acceptable.
  • min_max_rop_qty floor change can raise the minimum for SKUs whose sl_qty + fcst_lt − eoq/2 was previously below the old floor — verify no SKU is forced above its target-driven ROP by the floor in the A/B.
  • Transfer routes (eoq=0): buffers must be unchanged vs the old build (Type-1 fallback parity) — assert in the A/B.

Key files

File Role
files/MEIO/intial/Distributions.py V1 legacy Python (averaging loop, −eoq/2, EOQ-blind Normal inverse)
files/MEIO/intial/MEIO.py V1 legacy greedy optimizer (initialJump, further_jump, minMaxROPQty)
files/MEIO_v3/src/distributions.rs V2/V3 Rust forward fill_rate_for_rop + inverse *_fr_to_rop
files/MEIO_v3/src/marginal.rs V2/V3 Rust initial_jump, further_jump, min_max_rop_qty
files/MEIO_v3/src/math.rs eoq_step (V2 lot-size fix); target for the new loss module
files/inventory_engine/rop.py V3 Python reference (pure Type-1 percentile today)
files/meio_runner.py Loads eoq / lot_size / avg_size per SKU (:943, :1163) and dispatches the optimizer (:375)
.kilo/plans/1782658808592-meio-eoq-rop-loss-function.md Implementation plan (this chapter's companion)