Fun math algorithms

[Han]

Here are a few fun math tricks that I collected and still remember over the years. I'm curious what non-calculator math routines any of you may have learned (with or without explanation). I remember one more (computing square roots by hand which I did not learn until after grad school) which I will post if someone doesn't already do so after a few days.

Dividing by 9

Suppose we wish to compute 3461 / 9 by hand, without a calculator. The first digit of the answer is the first digit of the numerator: 3. The next digit in the answer comes from taking the previous digit, namely 3, and adding the next digit in the numerator modulo 9: \( 3+4 \pmod 9 \equiv 7 \pmod 9\). If the sum was greater than 9, add 1 to the previous digit. Since 7 is clearly less than 9, we do not add 1 to the preceding digit 3, and our partial answer is 37. Repeat the process of taking the current digit and adding the next digit in the numerator. So \( 7+6 \pmod 9 \equiv 13 \pmod 9 \equiv 4 \pmod 9\). Now the adjusted partial answer is 384 -- note that the second digit was previously 7, but since the third digit sum was 13 (greater than 9), we adjust the preceding digit by 1. Once we reach the one's digit in the numerator (in this case 1), the remainder would be 4+1 (4 from the previous digit and 1 from the next digit in the numerator). (If this sum were larger than 9, we would again increment the preceding digit by 1.) So the result is 384 with a remainder of 5.

To see how it works, it is best to study the following expansion by 10 and 9:
\[
\begin{align}
3461 & = 3\cdot 10^2 \cdot 10 + 4 \cdot 10 \cdot 10 + 6 \cdot 10 + 1 \\
& = 3\cdot 10^2 \cdot (9+1) + 4\cdot 10 \cdot 10 + 6\cdot 10 + 1 \\
& = 3\cdot 10^2 \cdot 9 + 3\cdot 10^2 + 4\cdot 10 \cdot 10 + 6\cdot 10 + 1\\
& = 3\cdot 10^2 \cdot 9 + (3+4)\cdot 10\cdot 10 + 6\cdot 10 + 1\\
& = 3\cdot 10^2 \cdot 9 + (3+4) \cdot 10 \cdot (9+1) + 6\cdot 10 + 1\\
& = 3\cdot 10^2 \cdot 9 + (3+4) \cdot 10 \cdot 9 + 7\cdot 10 + 6\cdot 10 + 1\\
& = 3 \cdot 10^2 \cdot 9 + (3+4)\cdot 10 \cdot 9 + (7+6)\cdot 10 + 1\\
& = 3\cdot 10^2 \cdot 9 + (3+4) \cdot 10 \cdot 9 + (9 + 4)\cdot 10 + 1\\
& = 3\cdot 10^2 \cdot 9 + (3+4) \cdot 10 \cdot 9 + 9\cdot 10 + 4\cdot 10 + 1\\
& = 3\cdot 10^2 \cdot 9 + (3+4+1) \cdot 10 \cdot 9 + 4\cdot 10 + 1\\
& = 3\cdot 10^2 \cdot 9 + (3+4+1) \cdot 10 \cdot 9 + 4\cdot (9+1) + 1\\
& = 3\cdot 10^2 \cdot 9 + (3+4+1) \cdot 10 \cdot 9 + 4\cdot 9 + (4 + 1)\\
\end{align}
\]
So
\[ \frac{3461}{9} = \frac{3\cdot 10^2 \cdot 9 + (3+4+1) \cdot 10 \cdot 9 + 4\cdot 9 + (4 + 1)}{9} =
3\cdot 10^2 + (3+4+1) \cdot 10 + 4+ \frac{4 + 1}{9}
\]

Squaring a number ending in 5:

Take the number formed by removing the final digit 5, multiply by 1 more than that number, and append 25 to the product. For example, \(6^2 = 4225\) where the first two digits 42 are obtained by \( 6\cdot (6+1) = 42\).

Explanation: write any number \(a\) ending in 5 as \( x5\). Then
\[ \begin{align}
a^2 & = (10\cdot x + 5)^2 \\
& = 100x^2 + 2\cdot 5 \cdot 10 \cdot x + 5^2\\
& = 100 (x^2+x) +25^2\\
& = 100\cdot x \cdot (x+1) + 25
\end{align} \]

Divisibility tests (I'll show divisibility by 11; you are welcome to try variations of this for other small divisors)

To test if a number n is divisible by 11, take alternating sums of the digits. If the sum is equivalent to 0 modulo 11, then the original number is divisible by 11. For example, consider 1254. Since \( -1 + 2 - 5 + 4 \pmod 11 \equiv 0 \pmod 11 \), the number 1254 is divisible by 11.

Note that \( 10^n \pmod{11} \equiv (-1)^n \pmod{11} \). Thus a number with digits \( d_9d_8d_7\cdots d_1 d_0 \) can be written as
\[ \begin{align}
d_9 \cdot 10^9 + d_8 \cdot 10^8 + \dotsm + d_1\cdot 10 + d_0 \pmod{11} & \equiv
d_9 \cdot (-1)^9 + d_8 \cdot (-1)^8 + \dotsm + d_1 \cdot (-1) + d_1 \pmod{11}\\
& \equiv -d_9 + d_8 - d_7 + \cdots + d_2 - d_1 + d_0 \pmod{11}
\end{align}
\]
This type of computation reminded me of a problem from a 7th grade mathematics competition where one had to determine whether a number (with many digits) written in base 13 was divisible by 12. Not only were calculators not allowed, but it was infeasible to expand the number into base 10 and then divide by 12 with the allotted time.

[telemachos]

From the Talmud, and respected as a rule for divisibility by 7 (remainder correct) by John Horton Conway:


5780? (5780/7 = 825, remainder 5)

57 80

114 + 80 = 194

1 94

2 + 94 = 96

96/7 = 13, remainder 5

See here: <https://arxiv.org/pdf/1903.04903.pdf>

[Albert Chan]

(09-05-2020 10:31 PM)Han Wrote:  Dividing by 9
\[ \frac{3461}{9} = \frac{3\cdot 10^2 \cdot 9 + (3+4+1) \cdot 10 \cdot 9 + 4\cdot 9 + (4 + 1)}{9} =
3\cdot 10^2 + (3+4+1) \cdot 10 + 4+ \frac{4 + 1}{9}
\]

Great trick ! Doing it vertically, you can clearly see the advantage

3461 ÷ 9
 76       // 3+4=7
 13       // 7+6=13
  41      // 1+3=4
   5      // 4+1=5
---------
384 + 5/9

Another trick is to scale denominator to 1 - ε, 1/(1-ε) = 1 + ε + ε² + ...

\(\frac{3461}{9} = \frac{34610+3461}{99} = \frac{380.71}{1-0.01} = 380.71 + 3.8071 + \cdots = 384 + 5/9\)

\(\frac{1254}{11} = \frac{12540-1254}{99} = \frac{112.86}{1-0.01} = 112.86 + 1.1286 + \cdots = 114\)

---

Instead of divisibility test, we might as well get the remainder.

1001 = 7×11×13, we can do all factors modulo together.
1000 ≡ -1 (mod 1001), so group digits in sets of 3.

123,456,789 (mod 1001) ≡ 789 - 456 + 123 ≡ 456

123,456,789 (mod 7 ) ≡ 456 ≡ 36 ≡ 1
123,456,789 (mod 11) ≡ 456 ≡ 16 ≡ 5
123,456,789 (mod 13) ≡ 456 ≡ 66 ≡ 1

[Han]

(09-06-2020 12:46 AM)Albert Chan Wrote:  Another trick is to scale denominator to 1 - ε, 1/(1-ε) = 1 + ε + ε² + ...

\(\frac{3461}{9} = \frac{34610+3461}{99} = \frac{380.71}{1-0.01} = 380.71 + 3.8071 + \cdots = 384 + 5/9\)

\(\frac{1254}{11} = \frac{12540-1254}{99} = \frac{112.86}{1-0.01} = 112.86 + 1.1286 + \cdots = 114\)

That is very nice!

[Albert Chan]

(09-06-2020 12:46 AM)Albert Chan Wrote:  Another trick is to scale denominator to 1 - ε, 1/(1-ε) = 1 + ε + ε² + ...

We can reduce summing terms, from O(n), to O(ln(n))

1/(1-ε) = 1 + ε + ε² + ... = (1+ε) + ε²(1+ε) + ε⁴[(1+ε) + ε²(1+ε)] + ...

For binary system, scaling by powers-of-2 can be done with shift in exponent.

If denominator is between 2/3 and 4/3, |ε| ≤ 1/3 → ε^32 < 5.4e-16
Thus, maximum of 5 sums is needed to reach IEEE double precision.

Example: 123/456 = (123/512) / (456/512) = 0.240234375 / 0.890625

>>> x, eps = 0.240234375, 1 - 0.890625
>>> x += x*eps; eps *= eps       # x = 0.266510009765625
>>> x += x*eps; eps *= eps       # x = 0.26969823986291885
>>> x += x*eps; eps *= eps       # x = 0.26973683658086722
>>> x += x*eps; eps *= eps       # x = 0.26973684210526305
>>> x += x*eps; eps *= eps       # x = 0.26973684210526316 = 123/456

We have quadratic convergence, but results is not self-correcting.
In other words, all calculations must be done with full precisions.

For arbitrary precision calculation of reciprocals of n, we can use Newton's method.
Note: we do not let f = 1/n - x, f' = -1, since the whole point is to avoid division.

f = n - 1/x       → f' = 1/x²       → x - f/f' = x + x*(1-n*x)

For 0.5 ≤ n < 1.0, guess x = 2.9-n-n work well, converged within 4 iterations.

Redo 123/456 as 0.240234375 * (1/0.890625)

>>> n = 0.890625
>>> x = 2.9 - n - n         # x = 1.1187499999999999
>>> x += x*(1-n*x)       # x = 1.1227923583984374
>>> x += x*(1-n*x)       # x = 1.1228070173524727
>>> x += x*(1-n*x)       # x = 1.1228070175438596 = 1/n
>>> 0.240234375 * x     # = 123/456
0.26973684210526316

[David Hayden]

Subtraction without ever having to subtract from a 2-digit number. I actually came up with this in 5th grade.

12345
-9876
------

For the first digit, the normal method is to borrow from the 10's column and compute (10+5)-6 to get 9. Instead subtract the lower number from the upper (6-5=1) and then subtract the result from 10: (10-1=9)

More generally, if you must borrow to subtract a-b, instead of computing 10+a-b, compute b-a, then 10-(b-a)

For my little head, I could do the 2 easy subtractions faster and more accurately than the one subtraction from a 2-digit number. Hmm, maybe I invented RISC computing for humans in 1973!

[Albert Chan]

How to estimate car payments ?

XCas> C := I*N / (1 - (1+I)^-N)       // C = |N*PMT/PV|, "compounding factor"
XCas> series(C,I,polynom)

\(1
+\frac{I(N+1)}{2}
+\frac{I^2 (N^2-1)}{12}
+\frac{I^3 (-N^2+1)}{24}
+\frac{I^4 (-N^4+20N^2-19)}{720}
+\frac{I^5 (N^4-10N^2+9)}{480}\)

For car payments, N is usually not too big, we can use: \(\large C ≈ {(IN + 3)^2 + 3\over 12}\)

This estimate does not require converting annual interest rate to monthly rate.
All calculations could be done on a 4-banger

Example, using U.S. 2020 numbers

Quote:The average loan amount for a new car in the first quarter of 2020 was $33,631,
with an average interest rate of 3.6% for a 60-month loan.

I*N = (12*I) * (N/12) = 3.6% * 5 = 18%
C-1 ≈ (IN)*(6+IN)/12 = 18% * 6.18/12 = 9.27%

Fixed 2 mode:
33631 Enter 60 /    → 560.52 (monthly payment, if no interest)
9.27 %              →  51.96 (monthly finance charge, estimated)
+                   → 612.48 (monthly payment, estimated)

For reference, car payments (exact) = $613.31

[EdS2]

That's rather good, but did I fail to follow something? I don't understand where the cubic comes from - please elaborate, if you can!

[Albert Chan]

Hi, EdS2

I was cutting off taylor series of C-1, eliminated O(I³) terms, then eliminated I only terms.

C-1 = I*(N+1)/2 + I²*(N²-1)/12 + ... ≈ (IN)/2 + (IN)²/12 = IN*(6+IN)/12

Keeping I only terms make better estimate, but more computations.
Update: to reduce computations, we ignore I² only term.

C-1 ≈ I*(N+1)/2 + I²*(N²-1)/12 = IN*(6+IN)/12 + I*(6-I)/12 ≈ IN*(6+IN)/12 + I/2

Redo example, PV=33631, I=3.6%/12=0.3% , N=5*12=60

C-1 ≈ IN*(6+IN)/12 + I/2 = 9.27% + 0.15% = 9.42%

Fixed 2 mode:
33631 Enter 60 /    → 560.52 (monthly payment, if no interest)
9.42 %              →  52.80 (monthly finance charge, estimated)
+                   → 613.32 (monthly payment, estimated)

Car payment over-estimated only 1 penny

[Albert Chan]

(10-17-2020 08:51 AM)EdS2 Wrote:  I don't understand where the cubic comes from - please elaborate, if you can!

My guess is you were really asking where C comes from.

You can think of car payments like this:
Instead of paying the loan, put car payments into a bank, collecting interest (same rate).
At the end of terms, take out the money, and pay off the loan, all at once.

At the end of terms, your bank deposit have grown (last deposit have no interest):

-FV = PMT*(1+I)^(N-1) + PMT*(1+I)^(N-2) + ... + PMT

  |FV/PMT| = 1 + (1+I) + ... + (1+I)^(N-1)
(1+I)|FV/PMT| = (1+I) + ... + (1+I)^(N-1) + (1+I)^N

Subtract 2 expressions, we have: I |FV/PMT| = (1+I)^N - 1

Since we put off paying the loan, it has grown too, FV = PV*(1+I)^N
Both FV must match, we can simply substitute:

I |PV*(1+I)^N/PMT| = (1+I)^N - 1
I |PV/PMT| = 1 - (1+I)^-N

C = |N*PMT/PV| = I*N / (1 - (1+I)^-N)

Nice thing about C is, for small N, it is relatively flat, C-1 ≈ I*(N+1)/2
(if I=0, we have C-1=0, or no financing cost, as expeced)

Once we calculated C, we can ignore interest rate: N |PMT| = C |PV|

[EdS2]

Thanks Albert. I'm shocked to notice - eventually - that I'd misread the expression for C as a cubic, when it's clearly quadratic. All is much clearer now!

[Albert Chan]

(10-17-2020 11:27 AM)Albert Chan Wrote:  C-1 ≈ IN*(6+IN)/12 + I/2

Let's rename C as \(C_+\), to signal compounding effect in the forward direction.

Knowing FV=0, we have \(N |PMT_+| = C_+ |PV|\)

For reverse direction, knowing PV=0, we have \(N |PMT_-| = C_- |FV|\)

\(C_- = \Large{C_+ \over (1+I)^N} = \frac{\left(I N \over 1-(1+I)^-N\right)}{(1+I)^N} = {I(-N) \over 1-(1+I)^N}\)

So, run the clock backwards (N → -N), \(C_+ \text{ turns into } C_-\)
Naturally, C₊ estimate formula apply:

\(C_± - 1 ≈ \large{IN (IN\;±\;6) \over 12} + {I\over2}\)

---
Lets build a car leasing formula !

\(N·PMT = N(PMT_+ - PMT_-) \)

PMT+ is car payments if we *buy* the car (not gives it back).
PMT− is credit of car payments, if we do give it back (turning buying, back to leasing)

\(N·PMT_+ = CAP · C_+ ≈ CAP \left(1 + \large{IN (IN\;+\;6) \over 12} + {I\over2}\right)\)
\(N·PMT_- = RES · C_- ≈ RES \left(1 + \large{IN (IN\;-\;6) \over 12} + {I\over2}\right)\)

Substract the 2, and divide by N, we have

\(PMT = \large \left(1 + {I\over2} + {(IN)^2 \over 12} \right) ({CAP - RES \over N})
+ ({I\over2}) \normalsize (CAP+RES) \)

[Albert Chan]

(10-19-2020 08:51 PM)Albert Chan Wrote:  \(PMT = \large \left(1 + {I\over2} + {(IN)^2 \over 12} \right) ({CAP - RES \over N})
+ ({I\over2}) \normalsize (CAP+RES) \)

If we remove the correction factor, we have the "official car lease formula"

Quote:PMT = (CAP – RES) ÷ TRM + (CAP + RES) x MF, where MF = APR/24 = I/2

Let's try some numbers: CAP = 30000, RES = 15000, 4% APR, 36 months lease.

Official leasing formula:
PMT = 15000/36 + 45000*(4%/24) ≈ 416.67 + 75 = 491.67

Add corrections:
I/2 + (IN)²/12 = (4%/24) + (12%)²/12 = 0.2866667%
Corrections = 416.67 * 0.2866667% ≈ 1.19
PMT = 491.67 + 1.19 = 492.86

Using HP50G Time Value of Money App:
PV = 30000, FV = -15000
N=36, APR=4%, P/YR=12, mode=END → PMT = -492.859775103

Car payment of $492.86, matching our estimate

[Albert Chan]

We can use \(C_±\) to solve Time value of money problem, if I is small, and N not too big.

XCas> C(d) := 1 + N*I*(N*I+d*6)/12 + I/2;       // compounding factor, d is direction
XCas> eqn := C(1)*PV + C(-1)*FV + N*PMT;       // variables follow cash flow sign convention

eqn is linear if solving for PV, FV, PMT; quadratic if solving for I,N

XCas> subst(eqn, [N, I, PV, FV] = [36, 0.04/12, 30000, -15000])
      → 36*PMT + 17743.0

For previous post example, car payments = $17743/36 = $492.86 (all digits correct)

What if car payments is $550, what is the APR ?

XCas> proot(subst(eqn, [N, PV, PMT, FV]=[36, 30000, -550, -15000])) * 12
      → [-6.12521299716, 0.0696574416048]

eqn estimated APR of 6.966% (all digits correct)