JavaScript implementation of Bernoulli numbers algorithm as described by Ada Lovelace

.gitignore

node_modules/
*~

README.md

Ada's Bernoulli Number Function

Function to calculate Bernoulli numbers using the same algorightm Ada Lovelace used 1842 in (what is widely regarded as) the first ever computer program. Derivation of the formula and the Analytical Engine implementation appear in Note G of Lovelace's notes on Italian mathematician Luigi Menabrea's essay published 1842.

Usage

First clone this gist and install dependencies

git clone https://gist.github.com/terotil/3f83a473f372d31f55d5 ada-bernoulli
cd ada-bernoulli
npm install

All the Bernoulli numbers are yours to take!

var adaBernoulliNumber = require("../bernoulli");
adaBernoulliNumber(5);
// => 1/42

Development

Run tests

npm test

bernoulli.js

/**
 * Expose `adaBernoulliNumber()`
 */
exports = module.exports = adaBernoulliNumber;

var Ratio = require("lb-ratio");

/**
 * Calculate n:th Bernoulli number using the same algorightm Ada
 * Lovelace used 1842 in the first ever computer program.
 *
 * Examples:
 *
 *     adaBernoulliNumber(0);
 *     // => -0.5
 *
 * @param {Number} n
 * @return {Ratio}
 * @api public
 */
function adaBernoulliNumber(n) {
  var termArgument, termBnIndex, termCount, termIndex, termSum, termVal;
  switch (false) {
    case !(n <= 0):
      throw "Positive index expected, got " + n;
    case !((n % 2) === 0):
      // All even (in Ada's indexing) Bernoulli numbers are zero
      return Ratio();
    default:
      // count of non-zero terms in addition to A0 in Ada's formula (8.)
      termCount = (n - 1) / 2;

      // value for which the n:th polynomial term equals 1 in Ada's formula (8.)
      termArgument = termCount + 1;

      // value of polynomial term A0
      termSum = Ratio(-1, 2).multiply(Ratio(2 * termArgument - 1,
                                            2 * termArgument + 1));

      // sum the polynomial term values of further terms
      for (termIndex = 0; termIndex <= termCount - 1; ++termIndex) {
        termBnIndex = 2 * termIndex + 1; // the index of term's Bernoulli number
        termVal = polynomialTerm(termBnIndex, termArgument);
        termSum = termSum.add(adaBernoulliNumber(termBnIndex).multiply(termVal));
      }
      // Finding Bn using the equation (8.) results in all the terms
      // on opposite side of Bn, thus the switch of sign.
      return termSum.negate().simplify();
  }
}

/**
 * Calculate value of the polynomial multiplier term of given `degree`
 * for given value of `n`.
 */
function polynomialTerm(degree, n) {
  var denominator, i, numerator;
  numerator = 1;
  denominator = 1;
  for (i = 0; i < degree; ++i) {
    numerator *= 2 * n - i;
    denominator *= 2 + i;
  }
  return Ratio(numerator, denominator).simplify();
}

package.json

{
  "name": "ada-bernoulli",
  "description": "Calculate Bernoulli numbers using the same algorightm Ada Lovelace used 1842 in the first ever computer program",
  "version": "0.0.1",
  "main": "./bernoulli.js",
  "private": true,
  "author": "Tero Tilus",
  "homepage": "https://gist.github.com/terotil/3f83a473f372d31f55d5",
  "license": "MIT",
  "scripts": {
    "test": "mocha",
    "lint": "jshint --show-non-errors *.js"
  },
  "dependencies": {
    "lb-ratio": "*"
  },
  "devDependencies": {
    "mocha": "*",
    "jshint": "*"
  }
}

test.js

var assert = require("assert");
var Ratio = require("lb-ratio");
var adaBernoulliNumber = require("./bernoulli");

// Copied from http://oeis.org/A027641 and http://oeis.org/A027642
var bernoulliNumbers = {
  0:  Ratio(      1,    1),
  1:  Ratio(     -1,    2),
  2:  Ratio(      1,    6),
  3:  Ratio(      0,    1),
  4:  Ratio(     -1,   30),
  5:  Ratio(      0,    1),
  6:  Ratio(      1,   42),
  7:  Ratio(      0,    1),
  8:  Ratio(     -1,   30),
  9:  Ratio(      0,    1),
  10: Ratio(      5,   66),
  11: Ratio(      0,    1),
  12: Ratio(   -691, 2730),
  13: Ratio(      0,    1),
  14: Ratio(      7,    6),
  15: Ratio(      0,    1),
  16: Ratio(  -3617,  510),
  17: Ratio(      0,    1),
  18: Ratio(  43867,  798),
  19: Ratio(      0,    1),
  20: Ratio(-174611,  330)
};


describe('adaBernoulliNumber', function(){
  var n;

  var testFor = function(n) {
    var abn, bn;
    abn = adaBernoulliNumber(n - 1); // Ada's indexing is off-by-one
                                     // compared to modern indexing of
                                     // Bernoulli numbers.
    bn = bernoulliNumbers[n];
    it('should return '+bn.toString()+' as '+n+'th Bernoulli number', function(){
      assert(abn.equals(bn), "expected "+bn.toString()+" got "+abn.toString());
    });
  };

  for (n = 2; n <= 20; ++n) {
    testFor(n);
  }
});

This is incredible! Great work!