Coffee Math

Feb 15, 2025

Introduction

I was drinking my normal cup (or three) of coffee the other morning and I reflected on some of the interesting math behind it all, as one does. As I thought, I realized that there’s some more interesting math behind how caffeine moves in our body than meets the eye, and that it would make an interesting topic for me to write about. So, here goes!

A Simplified Model of Caffeine, er… Elimination

When we consume coffee, or any caffeinated beverage, it stays in our body for a while and then leaves in, shall we say, a liquid form. Interestingly, caffeine doesn’t leave our body at a constant rate, but rather has a “half life”. This half-life can vary, but this source says 5 hours. That means that if there are CC mg of caffeine in your body at 9:00 am, at 2:00 pm there will only be C2\frac{C}{2} mg remaining. We can state this as a differential equation:

dCdt=kC\frac{dC}{dt} = kC

We can solve for C using some Calc 1 hand-waving:

dCC=kdt\frac{dC}{C} = kdt dCC=kdt=lnC=kt+A\int \frac{dC}{C} = \int kdt = \ln{C} = kt + A C=ekt+A=αektC = e^{kt + A} = \alpha e^{kt}

At t=0t=0, C(t)=αC(t) = \alpha, so we’ll rewrite it like this for clarity:

C(t)=C(0)ektC(t) = C(0) e^{kt}

Given some half-life hh, we know that from some starting point (say t=0t=0) that hh units of time later half the caffeine will remain in our body:

C(h)=C(0)2=C(0)ekhC(h) = \frac{C(0)}{2} = C(0) e^{kh} 12=ekh\frac{1}{2} = e^{kh} ln12=kh\ln{\frac{1}{2}} = kh ln2=kh-\ln{2} = kh k=ln2hk = \frac{-\ln{2}}{h}

Modeling Coffee Consumption

Accordingly to Mayo Clinic a cup of coffee has about 96mg of caffeine. I typically drink 3 cups of coffee a day (288mg), and space that over the course of 2(-ish) hours. That comes out to 144mghr\frac{144 mg}{hr} on average.

A Brief Heaviside Aside

In order to model coffee consumption, we’re going to need a new tool. This tool is the Heaviside step function which is defined as a piecewise function:

H(t)={0t<01t0H(t) = \begin{cases} 0 & t \lt 0 \\ 1 & t \geq 0 \\ \end{cases}

Here’s what it looks like graphed:

In a way, you can see this function as modeling a “switch turning on” at t=0t=0.

We can combine two of these functions to create a “impulse” function. That is, a function that “turns on” for a period τ\tau then “turns off”:

R(t;τ)=H(t)H(tτ)R(t; \tau) = H(t) - H(t - \tau)

Here’s an example for τ=2\tau = 2.

A Heaviside Drinking Function

Using this impulse function RR, we can scale it such that it models my caffeine intake. If I drink for 2 hours, τ=2\tau = 2. And, if we multiple RR by the amount of coffee I consume an hour we have a function that for any tt models the rate caffeine enters my system. We’ll call that function D(t)D(t) (“D” for “drink” of course!):

D(t)=144R(t;2)D(t) = 144 R(t; 2)

But we can generalize this a little better. If we drink nn cups of of some caffeined beverage with mm mg of caffeine over a period of τ\tau hours we end up with a caffeine ingestion rate of mnτ\frac{mn}{\tau} mg per hour. This gives

D(t)=mnτR(t;τ)D(t) = \frac{mn}{\tau} R(t; \tau)

or

D(t)=mnτ(H(t)H(tτ))D(t) = \frac{mn}{\tau} \left( H(t) - H(t - \tau) \right)

A New Differential Equation

We can now update our original differential equation to include caffeine intake along with our existing model for caffeine half-life:

dCdt=kC+mnτ(H(t)H(tτ))\frac{dC}{dt} = kC + \frac{mn}{\tau} \left( H(t) - H(t - \tau) \right)

A Numerical Solution

There’s a few methods out there, but we can most easily get a quick approximation of the behavior of that differential equation with Euler’s Method:

C(t+δt)=C(t)+δtdC(t)dtC(t + \delta t) = C(t) + \delta t \frac{dC(t)}{dt}

Here’s an interactive graph of that over a 24 hour period. Note that t=0t=0 corresponds to “when you first start drinking coffee”. So, a 7:00am cup of coffee would correspond to a bedtime of around t=16t = 16 for 8 hours of sleep.




So this is super fun to play with IMO, and it shows some cool properties of this system:

  • If you drink your coffee very fast, you get a bigger “jolt” because you give the “first” caffeine less time to leave your body before the last bit comes in.
  • If you drink your coffee slowly, you end up going to bed with more caffeine in your system because the last bit of caffeine has less time to leave your body before bedtime.

Some Final Thoughts and Potential Next Steps

Some remaining things to explore with this model would be to consider an algebraic solution to the differential equation. Laplace transforms are an ideal candidate for this sort of problem, and can potentially empower us to consider other models for how caffeine enters our body. For example, our current model is setup to handle a single day, but I typically drink coffee every morning. Could we model D(t)D(t) as periodic function? And if D(t)D(t) is periodic, could some sort of Fourier analysis help here?

And, as interesting as I think this model is, it only models caffeine in your body and not necessarily its effects. To model that better we’d have to consider some pharmalogical effects like:

  • Onset of Action which would potentially delay when the caffeine actually hits your bloodstream, not instantaneously upon drinking like we did here.
  • Duration of Action which would give us a better model of how the effect of caffeine changes over time in your body. There may be models for that as well.

There’s potentially other pharmacological things to consider that I’m unaware of, but those two factors alone could convert this into a reasonably complex system of equations. These systems may not be solvable by hand, and we’d have to resort to a numerical approach like Euler’s method or Runge-Kutta.

These all sound like topics for later writing… Until then!