Software 13.5.1 (EOQ and discounts). Example 13.1.2: the EOQ, orders a year and annual cost; the cost of three other order sizes relative to the optimum; then the all-units discount schedule, with the EOQ at each price, the best feasible order in each range, and its total annual cost including purchases:

>>> import numpy as np
>>> D, K, i = 12000, 500.0, 0.20  # demand a year, cost per order, rate
>>> def cost(Q, price):               # ordering + holding + purchase, a year
...     return K * D / Q + i * price * Q / 2 + price * D
...
>>> eoq = lambda price: np.sqrt(2 * K * D / (i * price))
>>> Q = eoq(12.0)
>>> round(Q, 1), round(D / Q, 2), round(cost(Q, 12) - 12 * D, 1)
(2236.1, 5.37, 5366.6)
>>> for q in (1500, 3000, 4472):           # cost relative to the optimum
...     print(q, round((cost(q, 12) - 12*D) / (cost(Q, 12) - 12*D), 4))
...
1500 1.0808
3000 1.0435
4472 1.25
>>> breaks = [(0, 12.0), (2000, 11.5), (5000, 11.0)]   # all-units discounts
>>> for lo, price in breaks:
...     q = max(eoq(price), lo)            # best order in this price range
...     print(price, round(eoq(price), 1), round(q, 1),
...           round(cost(q, price), 1))
...
12.0 2236.1 2236.1 149366.6
11.5 2284.2 2284.2 143253.6
11.0 2335.5 5000 138700.0

Reading the output. \(Q^* = 2236\), \(5.37\) orders a year, Rs. 5367. Orders of \(1500\) and \(3000\) cost \(8.1\%\) and \(4.4\%\) more, and \(4472 = 2Q^*\) costs \(25\%\) more, as the ratio formula predicts. With discounts, the EOQ at Rs. 11 is \(2336\), below the \(5000\) needed for that price, so the candidate is \(5000\); it costs Rs. 138 700 against Rs. 143 254 for the best order at Rs. 11.5, and is optimal.

MATLAB/Octave: sqrt(2*K*D./(i*price)) for the three EOQs.

Software 13.5.2 (The newsvendor). The critical fractile and order quantity for normal demand, checked by simulating \(400\,000\) days and searching a grid of orders:

>>> from scipy.stats import norm
>>> mu, sd = 200, 40                       # newsvendor: demand N(200, 40^2)
>>> price, c, salvage = 10.0, 6.0, 1.0
>>> cu, co = price - c, c - salvage        # underage and overage costs
>>> ratio = cu / (cu + co); Q = norm.ppf(ratio, mu, sd)
>>> round(ratio, 4), round(Q, 2)
(0.4444, 194.41)
>>> rng = np.random.default_rng(0); d = rng.normal(mu, sd, 400000)
>>> def profit(q):
...     sold = np.minimum(d, q)
...     return np.mean(price * sold + salvage * (q - sold) - c * q)
...
>>> qs = np.arange(150, 251, 5)
>>> best = qs[np.argmax([profit(q) for q in qs])]
>>> best, round(profit(Q), 2), round(profit(mu), 2)
(195, 657.54, 656.14)

Reading the output. The fractile \(0.444\) gives \(Q^* = 194.4\). The simulation’s best order on the grid of steps of \(5\) is \(195\), in agreement, and the expected profit at \(Q^*\) is \(657.5\), against \(656.1\) for ordering the mean demand \(200\): the profit curve is flat, like the EOQ cost, but its peak is below the mean because leftovers cost more than shortages here.

MATLAB/Octave: norminv(cu/(cu+co), mu, sd).

Software 13.5.3 (An \((s, S)\) policy by simulation). Periodic review with Poisson demand of mean \(10\) a period, ordering cost \(50\), holding cost \(1\) and backorder cost \(8\) per unit per period: the average cost of each \((s, S)\) on a grid over \(50\,000\) simulated periods (the same demands for every policy), the best pair, and the EOQ-based policy with \(S - s = \sqrt {2KD/h}\):

>>> rng = np.random.default_rng(1)
>>> dem = rng.poisson(10, 50000)           # demand per period
>>> Kf, h, p = 50.0, 1.0, 8.0              # order, holding, backorder costs
>>> def avg_cost(s, S):                    # periodic review, zero lead time
...     x, total = S, 0.0
...     for dt in dem:
...         if x <= s: total += Kf; x = S    # order up to S
...         x -= dt
...         total += h * max(x, 0) + p * max(-x, 0)
...     return total / len(dem)
...
>>> grid = [(s, S) for s in range(0, 16) for S in range(s + 1, 61, 2)]
>>> costs = {sS: avg_cost(*sS) for sS in grid}
>>> best = min(costs, key=costs.get); best, round(costs[best], 3)
((6, 35), 30.852)
>>> q = round(np.sqrt(2 * Kf * 10 / h)); q, round(avg_cost(5, 5 + q), 3)
(32, 30.917)

Reading the output. The best policy on the grid is \((6, 35)\), at an average cost of \(30.85\) a period: reorder when the stock falls to \(6\), which with zero lead time balances backorders against holding, and order up to \(35\), so that each order is at least \(35 - 6 = 29\) units, a little more on average because the stock usually falls below \(6\) before the review. The EOQ heuristic, \(S - s = 32\) with \(s = 5\), costs \(30.92\), within \(0.2\%\). Using common random numbers makes the comparison between policies far more precise than the absolute costs.

MATLAB/Octave: poissrnd(10, 1, 50000) and the same loop.