junosuarez · January 9, 2020 06:11 · Jan 9, 2020 · Jan 9, 2020
diff --git a/π.js b/π.js
@@ -72,12 +72,8 @@ function simulate(trials = 10) {
     if (result) {
       pass++;
     }
-    // process.stdout.moveCursor(50, -1);
-    // console.log(result);
   }
 
-  // console.log("-".repeat(40));
-
   /*
     Buffon's needle formula is p = 2L/Dπ
     We simulated p, so we can solve for π:

diff --git a/π.js b/π.js
@@ -0,0 +1,119 @@
+/*
+Approximating π by simulation 
+see live on runkit at https://runkit.com/junosuarez/5e16c3e9ccfea2001ae2c8f5
+
+In 1733 the French naturalist Georges Louis Leclerc, Comte de Buffon, posed and
+solved the following problem in geometric probability: when a needle of length L is
+dropped onto a wooden floor constructed from boards of width D (where D ≥ L),
+what is the probability that it will lie across a seam between two boards?
+Buffon determined that the probability is 2L/Dπ.
+via https://www.maa.org/sites/default/files/s1000042.pdf
+
+Can try to simulate this geometrically by trials of a random drop angle and origin and checking for
+intersections of board lines
+*/
+
+function trial(L, D) {
+  // boards are parallel to the y axis. so it becomes a vertical line test.
+  // the needle is a line from A to B with length L
+
+  // assert precondition
+  if (!(D >= L)) {
+    throw new Error("D must be >= L");
+  }
+
+  // let our scale be +- 10 board widths
+  const scale = 10 * D;
+
+  // first we generate a random point A and a random slope m
+  const A = [Math.random(), Math.random()].map(n => n * scale);
+  const Ax = A[0];
+  const m = Math.random();
+
+  // then we solve the triangle to find B from L
+  // we know the hypotenuse and the slope, need to calculate opp and adj sides
+  // we can shortcut because we only need to know Bx, not B
+  const theta = 2 * Math.PI * m;
+  // by law of cosines, ratio between hypotenuse and adjacent side is cosine theta
+  const Δx = Math.cos(theta) * L;
+  const Bx = Ax + Δx;
+
+  // by Intermediate Value Theorem, we pass the vertical line test if  Ax < Dn < Bx for any n in N
+  // console.log(Ax, "\t", Bx);
+
+  if (Ax === Bx) {
+    // if it's a vertical line, only pass if Ax = Dn for some n in N
+    return Ax % D === 0;
+  }
+
+  // naive iterate board widths, can be improved
+  // it doesn't matter which end of the line is A or B as long as the length is L
+  const min = Math.min(Ax, Bx);
+  const max = Math.max(Ax, Bx);
+  for (let x = 0; x < max; x += D) {
+    if (min < x && x < max) {
+      // console.log(min, "\t", x, "\t", max, "\t", true);
+      return true;
+    }
+  }
+
+  // console.log(min, "\t", "Dn", "\t", max, "\t", false);
+  return false;
+}
+
+function simulate(trials = 10) {
+  const L = 5;
+  const D = 10;
+
+  let pass = 0;
+
+  for (let i = 0; i < trials; i++) {
+    const result = trial(L, D);
+    if (result) {
+      pass++;
+    }
+    // process.stdout.moveCursor(50, -1);
+    // console.log(result);
+  }
+
+  // console.log("-".repeat(40));
+
+  /*
+    Buffon's needle formula is p = 2L/Dπ
+    We simulated p, so we can solve for π:
+    π = 2L/Dp
+  */
+  const p = pass / trials;
+  const π = (2 * L) / (D * p);
+
+  return { trials, p, π, e: Math.PI - π };
+}
+
+const results = [];
+for (let j = 0; j < 3; j++) {
+  let result = simulate(1e5);
+  results.push(
+    Object.assign({ run: j, ["abs e"]: Math.abs(result.e) }, result)
+  );
+}
+
+// display a pretty table of results
+const Table = require('easy-table')
+const t = new Table();
+results.forEach(r => {
+  Object.entries(r).forEach(([key, val]) => {
+    t.cell(key, val);
+  });
+  t.newRow();
+});
+t.sort(["abs e"]);
+t.total("abs e", {
+  printer: function(val, width) {
+    return "Avg: " + val;
+  },
+  reduce: function(acc, val, idx, len) {
+    acc = acc + val;
+    return idx + 1 == len ? acc / len : acc;
+  }
+});
+console.log(t.toString());