A faster way to calculate logarithm on the soroban

We will discuss a fast way to calculate logarithms on the Japanese abacus (soroban) up to the accuracy of 4 decimal places, without using logarithm tables. The method used here to produce an efficient algorithm for logarithm calculation can be applied to many other complicated functions.

Summa summarum we will implement the following formula, which has a maximum error lower than 10⁻⁵:

1 <= x <= 2

ln(x) ≈ -1.941142 + (3.529580 + (-2.461605 + (1.130888
	+ (-0.2888280 + 0.03111568 * x) * x) * x) * x) * x

This formula is very simple, using only addition, subtraction and multiplication. All the coefficients consist of 6 to 7 digits without the trailing zeros. They can even be memorized.

The final implementation on the soroban will have the following costs:

• Assignments: 6
• Additions: 1
• Subtractions: 5
• Multiplications: 6
• Divisions: 1

 

Table of contents

1. Method for producing efficient algorithms in Maple
2. Implementing the formula on the soroban
   1. Normalizing the value
   2. Calculation
   3. Full example

 

Method for producing efficient algorithms in Maple

There are not many resources available on the internet for logarithm calculation on the soroban. One of the best methods is detailed on this website: http://webhome.idirect.com/~totton/soroban/Logarithms/
It is still quite complicated and yields to the accuracy of only 3 digits.

After searching for a while I’ve found this idea from prof. Robert Israel: https://math.stackexchange.com/a/61347/403539

In the comments he gives a short explanation how that result was produced, which I was able to reproduce using Maple. The commands are the following:

with(numapprox):
Digits := 7
minimax(ln(x), x = 1..2, 5, 1, 'maxerror')
maxerror

Press enter after each line, so don’t copy-paste the whole text but line-by-line.
Click on the image for full size:

We didn’t get the same results, because a different precision (“Digits” directive) was chosen.
The important part is the “minimax” command which uses the Remez algorithm to find the best approximating polynomial for a function.

If we right-click on the expression and select “expand”, then we get the canonical form of the polynomial:

-1.941142 + 3.529580*x - 2.461605*x² + 1.130888*x³ - 0.2888280*x⁴ + 0.03111568*x⁵

Which makes it clear why we called the numbers coefficients. Although this form is not useful on the soroban, so we will use the expression with the parentheses.

The parameters of the “minimax” command are the following

1. “ln(x)”: The first parameter contains the expression to be approximated. It can be complicated, like: “sin(x² + 3*x)”, etc.
2. “x = 1..2″: The second parameter limits the approximation to the interval 1.0 <= x <= 2.0. If this interval is increased to 1..10, then the maxerror also increases, but around 100-fold.
3. “5″: The third parameter tells the maximum degree of the resulting polynomial. Which is also the maximum number of multiplications. Increase this if you want to reduce the error.
4. “1″: The fourth parameter is a weight function, with which one can ask for more precision for specific intervals and less for others. We used “1″ here, which asks for equal precision on all points of the 1..2 interval.
5. “‘maxerror’”: In the fifth parameter we can specify the name of a variable. The minimax command will return the maximal error in that variable.

Feel free to play with the parameters and with the other commands. Notice that increasing the “Digits” makes only a subtle difference in the accuracy. The maximal error can be decreased mostly by allowing polynomials of higher degrees. Also notice that Maple will return an error message if you try to generate a high-degree polynomial with a low “Digits” setting.

 

Implementing the formula on the soroban

Normalizing the value

Since the formula works only on the 1..2 interval, we will have to do something with larger numbers.

Basically we will always use the following decomposition. This formula makes use of the 2 most useful properties of the logarithm function: it converts multiplication to addition and the exponentiation to multiplication.

log(y) = log(2ⁿ * x) = n*log(2) + log(x)

Example value: 724.552 (We want to calculate the logarithm of this.)
Since it is outside of the 1..2 interval, we will have to normalize it.
Now we have to find the largest power of 2 which is still smaller (or equal) than our number:
Have a look at a list of powers of 2:
1, 2, 4, 8, 16, 32, 64, 128, 256, 512, 1024, 2048, 4096, 8192, …

As can be seen, our value is between these two powers: 512 <= 724.552 < 1024
So the power we were searching for is: 512 = 2⁹
We can get the x value by dividing 724.552 with 512:
724.552 / 512 = 1.415

log(724.552) = log(512 * 1.415) = log(2⁹ * 1.415) = 9*log(2) + log(1.415)

We will show how to calculate log(1.415) in the next paragraph. The value of log(2) can be memorized: ln(2) = 0.69315, log₁₀(2) = 0.30103

Note that there’s another way for the decomposition: by repeatedly dividing the value by 2 until we reach the 1..2 interval. In this case we have to remember (or store it on the other side of the soroban) the number of times we divided by 2. That number will be the power.

Calculation

This is the approximating formula already shown:

1 <= x <= 2

ln(x) ≈ -1.941142 + (3.529580 + (-2.461605 + (1.130888
	+ (-0.2888280 + 0.03111568 * x) * x) * x) * x) * x

The x value has to be normalized already in the range 1..2.

The raw algorithm of that formula is the following:

x should be already on the soroban as the result of the normalization. Then:

	1.	* 0.03111568
	2.	- 0.288828
	3.	* x
	4.	+ 1.130888
	5.	* x
	6.	- 2.461605
	7.	* x
	8.	+ 3.529580
	9.	* x
	10.	- 1.941142

But this can’t be used directly because of the negative values involved. Since using complementary values would be too complicated we will divide the soroban into two halves in the next chapter.

Full example

We want to calculate the natural logarithm of 724.552

Divide the soroban into 2 halves:
- The left side should store negative results
- The right side the positives

Always round the results to 6 decimal places after the comma.

	1.	Set up the value on the right side of the soroban:
			Right = 724.552
	2.	Search for the largest power of 2
		which is still smaller than 724.552
			2⁹ = 512
	3.	Right / 512 = 1.415 (write this number and 2⁹ down on a paper)

Now the value is normalized and we can start to evaluate the polynomial:

	4.	Right * 0.03111568 = 0.044029 (Coefficient of x⁵)
	5.	Left = 0.288828 (Coefficient of x⁴)
	6.	Left - Right = 0.244799
	7.	Left * 1.415 = 0.346391
	8.	Right = 1.130888 (Coefficient of x³)
	9.	Right - Left = 0.784497
	10.	Right * 1.415 = 1.110063
	11.	Left = 2.461605 (Coefficient of x²)
	12.	Left - Right = 1.351542
	13.	Left * 1.415 = 1.912432
	14.	Right = 3.529580 (Coefficient of x¹)
	15.	Right - Left = 1.617148
	16.	Right * 1.415 = 2.288265
	17.	Right - 1.941142 = 0.347123 (Coefficient of x⁰)

That is the result of ln(1.415). (accurate up to 4 decimal places)
Now we use the left side to store 9*ln(2):  (although it is a positive value)

	18.	Left = 0.693147  (this is ln(2))
	19.	Left * 9 = 6.238323	(because 2⁹ = 512)

To get the final result:

	20.	Right + Left = 6.585446

After rounding to 4 decimal places we get an accurate value:
	6.5854 = ln(724.552)

Notice that during the polynomial calculation we used only subtraction. That is because in the additions one of the operands was always a negative value.