Fix squareRoot returning overestimates (squareRoot(6) = 3)

src/algorithms/math/square-root/__test__/squareRoot.test.js

@@ -27,7 +27,7 @@ describe('squareRoot', () => {
     expect(squareRoot(0, tolerance)).toBe(0);
     expect(squareRoot(1, tolerance)).toBe(1);
     expect(squareRoot(2, tolerance)).toBe(1.4);
-    expect(squareRoot(3, tolerance)).toBe(1.8);
+    expect(squareRoot(3, tolerance)).toBe(1.7);
     expect(squareRoot(4, tolerance)).toBe(2);
     expect(squareRoot(15, tolerance)).toBe(3.9);
     expect(squareRoot(16, tolerance)).toBe(4);
@@ -66,4 +66,24 @@ describe('squareRoot', () => {
     expect(squareRoot(4.5, 10)).toBe(2.1213203436);
     expect(squareRoot(217.534, 10)).toBe(14.7490338667);
   });
+
+  it('should not overestimate the root', () => {
+    expect(squareRoot(6)).toBe(2);
+    expect(squareRoot(12)).toBe(3);
+    expect(squareRoot(35)).toBe(6);
+    expect(squareRoot(99)).toBe(10);
+    expect(squareRoot(6, 1)).toBe(2.4);
+    expect(squareRoot(12, 2)).toBe(3.46);
+  });
+
+  it('should match Math.sqrt for small, large and fractional numbers', () => {
+    expect(squareRoot(0.25, 10)).toBe(0.5);
+    expect(squareRoot(0.0001, 10)).toBe(0.01);
+    expect(squareRoot(1e12)).toBe(1e6);
+    expect(squareRoot(1e300)).toBe(1e150);
+
+    for (let n = 1; n <= 1000; n += 1) {
+      expect(squareRoot(n, 6)).toBe(Math.round(Math.sqrt(n) * 1e6) / 1e6);
+    }
+  });
 });

src/algorithms/math/square-root/squareRoot.js

@@ -17,22 +17,20 @@ export default function squareRoot(number, tolerance = 0) {
     return 0;
   }
 
-  // We will start approximation from value 1.
-  let root = 1;
+  // We will start approximation from value 1. After the first Newton step the
+  // approximation is never below the real root (AM-GM inequality), and every
+  // next step moves it down towards the root. So we iterate while it keeps
+  // decreasing, which gives the root with full floating point precision.
+  //
+  // Stopping when |number - root^2| <= delta is not enough: that bounds the error
+  // of the square, not of the root, and gives e.g. squareRoot(6) = 3.
+  let root = (1 + number) / 2;
 
-  // Delta is a desired distance between the number and the square of the root.
-  // - if tolerance=0 then delta=1
-  // - if tolerance=1 then delta=0.1
-  // - if tolerance=2 then delta=0.01
-  // - and so on...
-  const requiredDelta = 1 / (10 ** tolerance);
-
-  // Approximating the root value to the point when we get a desired precision.
-  while (Math.abs(number - (root ** 2)) > requiredDelta) {
-    // Newton's method reduces in this case to the so-called Babylonian method.
-    // These methods generally yield approximate results, but can be made arbitrarily
-    // precise by increasing the number of calculation steps.
-    root -= ((root ** 2) - number) / (2 * root);
+  // Newton's method reduces in this case to the so-called Babylonian method.
+  let nextRoot = (root + (number / root)) / 2;
+  while (nextRoot < root) {
+    root = nextRoot;
+    nextRoot = (root + (number / root)) / 2;
   }
 
   // Cut off undesired floating digits and return the root value.