Stack-based square and cube roots

To calculate square roots mentally, we can make an educated guess and do one iteration of the Babylonian method, in which we average our guess with the quotient of the original number and our guess. Expressed mathematically, where S is the number we need the square root of and x is our initial guess, our first-order approximation will be

To make our initial guess, we'll use an algorithm that's easier to work mentally than express in code:

1. Split the number into pairs of digits, starting from the right.
2. Replace each group after the leftmost group with a zero.
3. Find the closest perfect square to the leftmost group.
4. Replace it by its square root.

Making a guess for 52,389, then,

starts the algorithm at 200.

The numerical code to do the same is somewhat opaque:

(defn guess-sqrt [x]
  (loop [x (float x)
         i 0]
    (if (< x 100)
      (* (int (Math/pow 10 i))
         (Math/round (Math/sqrt x)))
      (recur
        (/ x 100)
        (inc i)))))

To use this function as a stack operation, we need a monadic complement to dyad from the introductory post in this series:

(defn monad [f]
  ^:stack-op
  (fn [stack]
    (conj (pop stack) (f (peek stack)))))

Then we can wrap guess-sqrt with

(alter-var-root #'guess-sqrt monad)

Now we can execute one iteration of the Babylonian method:

(execute [52389] [dup guess-sqrt dup [/] dip avg])

[230.9725]

Less than 1% error from the actual value of 228.9. Walking through individual steps with executions shows that none of the calculations is particularly onerous to do mentally:

(executions
  [52389]
  [dup guess-sqrt dup [/] dip avg])

([52389]
 [52389 52389]
 [52389 200]
 [52389 200 200]
 [52389 200 200 [#function[clojure.core//]]]
 [261.945 200]
 [230.9725])

Answers don't have to be exact. In practice I frequently round calculations, such as halving 520 instead of 523. When greater precision is needed, I've found it convenient to work in fractions rather than decimals.

To generate an algorithm for calculating cube roots, we need to understand that the Babylonian method is a specific application of Newton's method, in which we're trying to find the zeros of the function f(x) = x² – S:

Extending this approach to cube roots gives us the formula

To generate an initial guess, we'll group digits in threes, pick a perfect cube from the leftmost group, and take its cube root:

In code,

(defn guess-cube-root [x]
  (loop [x (float x)
         i 0]
    (if (< x 1000)
      (* (int (Math/pow 10 i))
         (Math/round (Math/pow x (/ 3))))
      (recur
        (/ x 1000)
        (inc i)))))

(alter-var-root #'guess-cube-root monad)

The stack algorithm for the equation derived above is

(execute [52389]
  [dup guess-cube-root dup [sqr /] dip double + third])

[37.581043]

And step by step:

(executions [52389]
  [dup guess-cube-root dup [sqr /] dip double + third])

([52389]
 [52389 52389]
 [52389 40]
 [52389 40 40]
 [52389
  40
  40
  [#function[chunking-stack-math-operations/sqr]
   #function[clojure.core//]]]
 [32.743126 40]
 [32.743126 80]
 [112.743126]
 [37.581043])

Hiding within the execution of dip is an imposing division of 52,389 by 1,600, but practiced mental mathematicians will see the opportunity to estimate it with 512 ÷ 16 = 32, just over 2% from the exact quotient of 32.743126.

Personally I've found these stack algorithms to be harder to remember than the traditional mathematical notation. Standard math notation creates a two-dimensional, highly visual thought picture, while stack algorithms are linearized and verbal. However, performing the calculations as stack algorithms does demonstrate that no great feats of working memory are necessary to find surprisingly close estimates for square and cube roots.