Eat Now or Store for Later?
Every village, colony, and game-economy designer eventually faces the same dilemma. A harvest has just come in. Do you feed everyone generously today, which makes the population grow, or do you lock part of the grain away for a future that might be hungrier? Storage is not free, and some of the grain will rot. In this article we turn that dilemma into a precise optimization problem, solve it in Python with an exact gradient (so it runs fast), and visualize the answer in a single dashboard that includes 3D plots.
1. The Problem
A settlement lives for $T = 20$ turns. In every turn:
- A harvest $H_t$ arrives. It varies with a seasonal cycle on top of a slowly growing trend.
- The settlement decides how much of its available food to eat now and how much to store.
- Eaten food determines population growth. Stored food loses a fraction to spoilage and costs a fraction in upkeep.
The objective is to maximize the total population summed over all 20 turns.
The tension is easy to state:
- Eating more now grows the population sooner, and a larger population counts in every later turn too.
- But a larger population also needs more food later.
- Eating too much in a bumper-harvest turn wastes food, because growth saturates.
- Starving in a lean turn is punished heavily, because growth collapses.
- Stored food is a buffer against lean turns, but it decays.
2. Mathematical Model
State and decision. In turn $t = 0, 1, \dots, T-1$ the state is the population $P_t$ and the food stock $S_t$. The decision is the eating fraction $u_t \in [0, 1]$.
Food balance.
$$
A_t = S_t + H_t, \qquad c_t = u_t A_t, \qquad r_t = \frac{c_t}{P_t}
$$
Here $A_t$ is the available food, $c_t$ is the food eaten, and $r_t$ is the per-capita ration, where $r_t = 1$ means exactly one unit of food per person.
Population growth. The growth factor depends on the ration through a smooth saturating curve:
$$
P_{t+1} = P_t , g(r_t), \qquad g(r) = g_{\mathrm{lo}} + \left(g_{\mathrm{hi}} - g_{\mathrm{lo}}\right)\left(1 - e^{-(r/r_c)^2}\right)
$$
With no food the population shrinks by the factor $g_{\mathrm{lo}}$. With plentiful food it grows by at most $g_{\mathrm{hi}}$. The squared exponent makes the curve S-shaped, so mild shortages are tolerable but severe shortages are disastrous, and extra food beyond the knee of the curve brings little benefit.
Storage with spoilage and cost. Each stored unit independently spoils with probability $\delta$ per turn. For a large stock, the law of large numbers makes the expected stock an excellent description, so spoilage removes the fraction $\delta$ of what is stored. Storage upkeep costs $\kappa$ food units per stored unit per turn. Writing $\rho = 1 - \delta - \kappa$ for the survival factor:
$$
S_{t+1} = \rho , (1 - u_t) , A_t
$$
Objective.
$$
\max_{u_0, \dots, u_{T-1} \in [0,1]} ; J = \sum_{t=0}^{T-1} P_{t+1}
$$
Food left in storage after the final turn is worthless, so the model automatically pushes the settlement to run its stock down toward the end.
3. Concrete Example
| Quantity | Symbol | Value |
|---|---|---|
| Number of turns | $T$ | 20 |
| Initial population | $P_0$ | 100 |
| Initial stock | $S_0$ | 40 |
| Growth factor with no food | $g_{\mathrm{lo}}$ | 0.70 |
| Maximum growth factor | $g_{\mathrm{hi}}$ | 1.12 |
| Ration scale | $r_c$ | 1.0 |
| Harvest | $H_t$ | $(110 + 5t)\left(1 + a \sin\frac{2\pi t}{5}\right)$ |
| Seasonal amplitude | $a$ | 0.5 |
| Spoilage probability | $\delta$ | 0.15 |
| Storage cost | $\kappa$ | 0.03 |
The harvest grows linearly and oscillates with a period of five turns, so there are rich turns and lean turns. We compare three policies:
- Optimal: the best eating fractions found by the optimizer.
- Eat-all: consume everything every turn ($u_t = 1$).
- Constant 80%: always eat 80% of the available food and store the rest.
We also repeat the optimization over a grid of spoilage rates $\delta$ and seasonal amplitudes $a$ to see when storage is worth it.
4. Solution Method and Speed-Up
The problem has 20 continuous decision variables with simple bounds, and it is not convex. The growth curve is S-shaped, so a single starting point may end in a poor local optimum. We therefore use bounded L-BFGS-B from several different starting policies and keep the best result.
The naive approach lets the optimizer estimate gradients by finite differences, which costs 21 full simulations per gradient. We instead compute the exact gradient with a backward (adjoint) recursion. Let $J_t = \sum_{k=t}^{T-1} P_{k+1}$ be the population total from turn $t$ onward, and define the sensitivities
$$
\mu_t = \frac{\partial J_t}{\partial P_t}, \qquad \nu_t = \frac{\partial J_t}{\partial S_t}, \qquad \mu_T = \nu_T = 0 .
$$
With $a_t = 1 + \mu_{t+1}$ and $g’(r) = (g_{\mathrm{hi}} - g_{\mathrm{lo}}), e^{-(r/r_c)^2}, \dfrac{2r}{r_c^2}$, the recursion runs from $t = T-1$ down to $t = 0$:
$$
\frac{\partial J}{\partial u_t} = a_t, g’(r_t), A_t - \nu_{t+1}, \rho, A_t
$$
$$
\mu_t = a_t \left( g(r_t) - r_t, g’(r_t) \right), \qquad
\nu_t = a_t, g’(r_t), u_t + \nu_{t+1}, \rho, (1 - u_t)
$$
One forward pass plus one backward pass gives the whole gradient. That is roughly twenty times cheaper than finite differences, which makes the 121-point parameter sweep in Section 7 finish quickly.
5. Full Source Code
1 | import math |
6. Code Walkthrough
Parameters and helper functions. The constants at the top mirror the table in Section 3. harvest(amp) builds the harvest sequence, a linear trend multiplied by a sine wave whose size is controlled by amp. growth(r) is the vectorized growth curve used for plotting. growth_with_slope(r) is the scalar version that returns both $g(r)$ and $g’(r)$, which the gradient needs. It uses math.exp, which is faster than NumPy for single numbers and quietly returns zero for very negative arguments, so extreme rations cannot cause an overflow.
simulate. This is a literal translation of the model. For each turn it computes the available food, the amount eaten, and the ration, then updates the population with the growth factor and the stock with the survival factor $\rho$. It also returns the rations and consumption so we can plot them later.
neg_value_and_grad. This is the heart of the speed-up. The first loop is the same simulation, but it stores the intermediate growth factors, their slopes, the available food, and the rations. The second loop walks backward from the last turn and carries the two sensitivities wP and wS, the code versions of $\mu$ and $\nu$. At each turn it writes the exact derivative of the total population with respect to that turn’s eating fraction. The function returns the negated value and gradient because SciPy minimizes.
solve. Because the growth curve is S-shaped, the problem can have several local optima. The function therefore runs L-BFGS-B from three starting policies (eat everything, eat 80%, eat half) and keeps the best result. All eating fractions are bounded to $[0, 1]$ through the bounds argument. Passing jac=True tells SciPy that our function returns the value and gradient together, so no finite differences are needed.
Base case. The script solves the 20-turn problem for the base parameters, then simulates the two baseline policies for comparison. It prints a turn-by-turn table (harvest, stock, food eaten, eating fraction, ration, and ending population) and the three totals.
Parameter sweep. Two nested loops solve the full optimization on an 11 × 11 grid of spoilage rates (0% to 50%) and seasonal amplitudes (0 to 0.8). For every grid point we record the optimal total and the eat-all total, then convert the difference into a percentage gain. The exact gradient is what keeps this loop fast: 121 optimizations, each with three starts, finish in a short time.
Dashboard. All six panels are drawn in a single figure so that one image captures the whole story. The 3D panels use projection="3d", and view_init sets the camera angle so that the surfaces are easy to read.
7. Results
Console Output
============================================================================== Optimal policy (spoilage 15%, storage cost 3%, seasonal amplitude 0.5) ============================================================================== Turn Harvest Stock Eaten Eat frac Ration Pop(end) 1 110.0 40.0 136.6 0.910 1.366 105.5 2 169.7 11.0 130.4 0.722 1.237 108.5 3 155.3 41.2 121.3 0.617 1.117 108.5 4 88.3 61.7 110.4 0.736 1.018 105.3 5 68.2 32.4 100.6 1.000 0.955 100.2 6 135.0 0.0 135.0 1.000 1.347 105.4 7 206.6 0.0 144.6 0.700 1.372 111.3 8 187.6 50.8 139.9 0.587 1.257 115.0 9 105.9 80.8 133.0 0.712 1.156 116.1 10 81.3 44.1 125.4 1.000 1.080 114.9 11 160.0 0.0 160.0 1.000 1.393 121.7 12 243.5 0.0 170.5 0.700 1.401 129.1 13 220.0 59.8 164.9 0.589 1.277 134.0 14 123.6 94.2 155.9 0.716 1.164 135.6 15 94.4 50.7 145.1 1.000 1.071 133.7 16 185.0 0.0 185.0 1.000 1.383 141.5 17 280.4 0.0 232.7 0.830 1.645 154.5 18 252.3 39.1 224.4 0.770 1.452 165.2 19 141.2 54.9 196.2 1.000 1.188 168.1 20 107.5 0.0 107.5 1.000 0.640 141.3 ------------------------------------------------------------------------------ Total population over 20 turns Optimal policy : 2515.4 Eat-all policy : 2329.4 Constant 80% : 2303.2 Gain of optimal over eat-all: 7.98 % ------------------------------------------------------------------------------ Sweep: largest gain over eat-all = 39.83 % (spoilage 0%, amplitude 0.8) Sweep: smallest gain over eat-all = 0.00 %
Result Image
8. How to Read the Results
Population panel (top left). The three curves show how the same village evolves under three different philosophies. The legend reports each policy’s 20-turn sum, which is exactly the quantity we maximized. The gap between the cyan optimal curve and the red eat-all curve is the value created purely by timing the consumption. Whatever the exact numbers, the optimal curve cannot lie below the eat-all curve in total, because eat-all is one of the starting points of the search.
Harvest, consumption, and stock panel (top middle). The blue bars are the harvest, the cyan line is the food actually eaten, and the yellow dashed line is the stock at the beginning of each turn. Where the cyan line is smoother than the bars, the settlement is turning a spiky harvest into steady meals by moving food through the warehouse. The stock line should rise after rich harvests and fall during lean turns. Near the final turns it should shrink toward zero, because leftover grain at the end of turn 20 earns nothing.
3D trajectory panel (top right). Each policy is a path through the space of turn, population, and stock. The eat-all path hugs the floor of the stock axis because it never stores anything. The optimal path lifts off the floor whenever it builds a buffer and returns to it as the horizon approaches. Rotating the view in a notebook makes it easy to see how the buffer and the population interact.
Optimal total population surface (bottom left). This surface shows the best achievable total as a function of the spoilage rate and the seasonal amplitude. Higher spoilage means a leakier warehouse, which lowers what an optimal planner can achieve. The red dot marks the base case from Section 3 so the single run can be located on the surface.
Value of storing surface (bottom middle). This is the most decision-relevant picture. It plots the percentage improvement of the optimal policy over eat-all. The surface answers the practical question of when storage is worth the trouble. Storage should matter most where seasons swing widely and spoilage is low, since that is where moving food in time pays the most. It should matter least where the harvest is steady or spoilage is severe. The sweep summary lines in the console print the largest and smallest gains on the grid.
Growth curve panel (bottom right). The white curve is the S-shaped growth factor $g(r)$, and the colored dots are the rations the optimal policy actually chooses, colored by turn. Dots clustered near the knee of the curve reveal the central idea of the solution: because the curve flattens at high rations and falls steeply at low rations, an optimal planner avoids both overeating and starvation and keeps the ration in the productive middle zone. Any dots that stray far from that zone typically appear at the very end of the horizon, where the planner spends the last of its stock.
9. Takeaways
The model is small, but it captures the core logic of storage decisions:
- Smoothing is valuable because the growth response is nonlinear. Because starvation hurts more than abundance helps, shifting food from rich turns to lean turns raises the total even when part of the stored food is lost.
- Spoilage and upkeep set the price of smoothing. Each unit moved through the warehouse is discounted by $\rho = 1 - \delta - \kappa$, so storage is worthwhile only when the benefit of smoothing exceeds that loss.
- The horizon matters. With a finite horizon, the stock has no value at the end, so the optimal policy ends by consuming everything.
- Exact gradients make the exploration cheap. The adjoint recursion turned a 121-point parameter sweep into a task that finishes quickly, which lets us study the whole landscape instead of a single case.
Change $\delta$, $\kappa$, the harvest cycle, or the shape of the growth curve and rerun the script. Each change tells a different story about when it pays to eat now and when it pays to wait.
