Every Lunation Ends Where It Began

code-a-dayalgorithmsspacecalendarsmathematics

To generate a year of lithophane lunar phases, I need to know what the moon looks like on any given date. Not approximately — precisely enough to pick the right photograph from my library of 60-odd lunar captures. The calculation is ancient, but the arithmetic underneath it is one of computing’s most fundamental operations: the remainder.

The synodic month — new moon to new moon — runs 29.53059 days. January 6th, 2000 at 18:14 UTC was a new moon, a reference point astronomers call an epoch. From there, every moon phase becomes a subtraction and a division:

EPOCH = Time.utc(2000, 1, 6, 18, 14)
SYNODIC = 29.53059

def moon_phase(date)
  days = (date - EPOCH) / 86400.0
  age = days % SYNODIC
  illumination = (1 - Math.cos(2 * Math::PI * age / SYNODIC)) / 2.0
  [age.round(2), (illumination * 100).round(1)]
end

puts moon_phase(Time.utc(2026, 2, 15))
# => [5.78, 17.5]

The modulo operator (%) does the work that matters. Without it, you’d need a lookup table stretching from 2000 to 2026 — nearly ten thousand entries. With it, any date collapses into a position between 0 and 29.53.

The C version makes the arithmetic explicit, which is instructive because floating-point modulo isn’t a single CPU instruction the way integer modulo is:

#include <stdio.h>
#include <math.h>

int main(void) {
    double epoch = 946836840;  // 2000-01-06T18:14:00 UTC
    double now = 1771113600;   // 2026-02-15T00:00:00 UTC
    double synodic = 29.53059 * 86400;
    double age = fmod(now - epoch, synodic) / 86400;
    double illum = (1 - cos(2 * 3.14159265 * age / 29.53059)) / 2;
    printf("Age: %.2f days, %.1f%% illuminated\n", age, illum * 100);
    return 0;
}

fmod() is the floating-point remainder — slower than integer %, but necessary when your divisor isn’t a whole number.

The elegance is that modular arithmetic handles time’s forward march by pretending it’s circular. The moon doesn’t care about Unix timestamps or Julian dates. It orbits. The same calculation that ran on a PDP-11 in 1975 runs on my laptop today, because the cycle hasn’t changed and the remainder operation hasn’t changed. My script spits out [5.78, 17.5] for February 15th — a waxing crescent, 17.5% illuminated, the kind of slender arc that prints beautifully thin at the edges and thick toward the terminator.

Now I just need to match that number to a photograph. Which means sorting through 60 images by measured phase angle. Another script, another loop, another modulo lurking somewhere in the comparison logic.

The moon keeps going. The numbers wrap.