Friday 10 December 2010

Computus 2

Of course, if you're as old as I am, you'll remember programming on machines with no native division instruction (neither integer nor floating-point), so all those divisions and modulos aren't going to get you very far on a "old school" 16-bit computer. Using the information from Hacker's Delight, you can re-cast these operations to multiplications and shifts:

void easter_u16_nodiv(u16 year, u16* month, u16* day)
{
    u16 a = year >> 1;
    u16 b = year >> 2;
    u16 c = year >> 4;
    u16 d = b + c;
    u16 e = (a + d + (d >> 4) + (d >> 5) + (d >> 10)) >> 4;
    u16 f = year - e * 19;
    u16 g = f - ((f + 13) >> 5) * 19;
    d = (b >> 1) + (b >> 3);
    e = (d + (d >> 6) + (d >> 7)) >> 4;
    f = b - e * 25 + 39;
    d = e + (f >> 5);
    u16 h = (d * 3) >> 2;
    d = (d * 8) + 5;
    e = (d >> 1) + (d >> 3);
    e = (e + (e >> 6) + (e >> 7)) >> 4;
    f = d - (e * 25) + 7;
    d = (g * 19) + h - e - (f >> 5) + 15;
    e = ((d >> 1) + (d >> 5)) >> 4;
    f = d - (e * 30);
    d = f - ((f + 2) >> 5) * 30;
    g = d + ((29578 - g - (d * 32)) >> 10);
    f = b - h + g + 2;
    d = a + c + (f >> 1) + (f >> 4);
    e = (d + (d >> 6) + (d >> 12)) >> 2;
    f = year + f - (e * 7);
    e = g - f + ((f * 2341) >> 14) * 7;
    d = e >> 5;
    *day = e - d * 31;
    *month = d + 3;
}

Again, if you want to verify the algorithm (only valid for Gregorian years 0 to 65535), you can do so with Dr J R Stockton's test harness with the following Javascript snippet:

    var a, b, c, d, e, f, g, h;
    a = YR >>> 1;
    b = YR >>> 2;
    c = YR >>> 4;
    d = b + c;
    e = (a + d + (d >>> 4) + (d >>> 5) + (d >>> 10)) >>> 4;
    f = YR - e * 19;
    g = f - ((f + 13) >>> 5) * 19;
    d = (b >>> 1) + (b >>> 3);
    e = (d + (d >>> 6) + (d >>> 7)) >>> 4;
    f = b - e * 25 + 39;
    d = e + (f >>> 5);
    h = (d * 3) >>> 2;
    d = (d * 8) + 5;
    e = (d >>> 1) + (d >>> 3);
    e = (e + (e >>> 6) + (e >>> 7)) >>> 4;
    f = d - (e * 25) + 7;
    d = (g * 19) + h - e - (f >>> 5) + 15;
    e = ((d >>> 1) + (d >>> 5)) >>> 4;
    f = d - (e * 30);
    d = f - ((f + 2) >>> 5) * 30;
    g = d + ((29578 - g - (d * 32)) >>> 10);
    f = b - h + g + 2;
    d = a + c + (f >>> 1) + (f >>> 4);
    e = (d + (d >>> 6) + (d >>> 12)) >>> 2;
    f += YR - (e * 7);
    e = g - f + ((f * 2341) >>> 14) * 7;
    return e;

No comments:

Post a Comment