Sunday 13 November 2016

arithmetic - How does this algorithm calculates the natural logarithm?



This is the pseudocode



(0)

x = input
LOG = 0
while x >= 1500000:
LOG = LOG + 405465
x = x * 2 / 3
(1)
x = x - 1000000
y = x
i = 1
while i < 10:

LOG = LOG + (y / i)
i = i + 1
y = y * x / 1000000
LOG = LOG - (y / i)
i = i + 1
y = y * x / 1000000
return(LOG)


The algorithm uses integer numbers and it assumes that every number is multiplied by $10^6$ so the last $6$ digits represent the decimal part. For example, the number $20.0$ would be represented as $20 \times 10^6$.




I have understood how part ($1$) of the pseudocode works. It just uses the following Taylor Series:
$$\ln(x) = (x-1) - \frac{(x-1)^2}{2} + \frac{(x-1)^3}{3}-\frac{(x-1)^4}{4} +\ldots,$$



which is valid for $|x-1| \le 1$.



I don't understand though how the first while at ($0$) works. I know that its purpose is to get an $x$ such that $|x-1| \le 1$ but how does it calculate the value LOG? The number $405465$ does not seem to be an arbitrary number, so how was it chosen? As for the value 2/3 could that be changed with a different value provided that I change 1500000 accordingly?


Answer




The number 405465 does not seem to be an arbitrary number, so how was it chosen?





The natural logarithm of $\frac32$ is approximately $0.40546511$. Round that to six digits, scale up by the $10^6$ factor, and there it is.
So what does this do? Well, consider what happens if we input $\frac32$ as the number to find the logarithm - which scales to $x=1500000$. We go through the first loop once and add $405465$ to LOG, while multiplying $x$ by $\frac23$ to get $1000000$. That's low enough to terminate the loop.
Now, in the second loop, we reduce $x$ by $1000000$, dropping it to zero. We run the power series, adding $0$ to LOG ten times.



So, then, that number we added in the first loop is all we get for the logarithm of $\frac32$. It has to be that logarithm for the result to be correct.




As for the value $\frac23$ could that be changed with a different value provided that I change 1500000 accordingly?





We have three values to change together; the factor $1/k$ to multiply by, the threshold $10^6\cdot k$ to exit the first loop, and the value $10^6\cdot\ln(k)$ to add each time through the first loop. There are of course tradeoffs involved in time cost and/or accuracy. For larger $k$, we need more terms in the power series to maintain accuracy, while for smaller $k$ we take more trips through the first loop.


No comments:

Post a Comment

real analysis - How to find $lim_{hrightarrow 0}frac{sin(ha)}{h}$

How to find $\lim_{h\rightarrow 0}\frac{\sin(ha)}{h}$ without lhopital rule? I know when I use lhopital I easy get $$ \lim_{h\rightarrow 0}...