An attempt for approximating the logarithm function $\ln(x)$: Could be extended for big numbers?

264 Views Asked by At

An attempt for approximating the logarithm function $\ln(x)$: Could be extended for big numbers?

PS: Thanks everyone for your comments and interesting answers showing how currently the logarithm function is numerically calculated, but so far nobody is answering the question I am asking for, which is related to the formula \eqref{Eq. 1}: Is it correctly calculated?, Could a formula for the logarithm of large numbers be found with it? Here with "big/large numbers" I am meaning in the same sense of how the Stirling's approximation formula approximates the factorial function at large values.


Intro__________

On a previous question I found that the following approximation could be used: $$\ln\left(1+e^x\right)\approx \frac{x}{1-e^{-\frac{x}{\ln(2)}}},\ (x\neq 0) \quad \Rightarrow \quad \dfrac{\ln\left(1+x^{\ln(2)}\right)}{\ln\left(x^{\ln(2)}\right)} \approx \frac{x}{x-1}$$

And later I noted that I could do the following: $$\dfrac{\ln\left(1+x^{\ln(2)}\right)}{\ln(2)} \approx \frac{x\ln\left(x\right)}{x-1}$$ by differentiating both sides: $$\frac{x^{\ln(2)-1}}{1+x^{\ln(2)}} \approx \frac{1}{x-1}-\frac{\ln(x)}{(x-1)^2}$$

From where I could isolate the logarithm as an approximation function: $$\ln(x) \approx f_1(x) = \frac{\left(x-1\right)\left(x+x^{\ln(2)}\right)}{x\left(1+x^{\ln(2)}\right)}$$

which end been a quite good approximation: sacrificing little accuracy it follows the logarithmic curve for longer than Taylor expansions or Padé approximants of similar order: as example the 2nd order Taylor series center at $x=1$ equals $T(x)=(x-1)-\frac12 (x-1)^2$, and the Pade approximant of 2nd order in both numerator and denominator at $x=1$ is $P(x)=\frac{3(x^2-1)}{x(x+4)+1}$ are shown in the following figures 1 and 2 compared with $f_1(x)$:

simple approximation

It follows the graph quite good and way longer that the classic approximations just with some little crossovers, but as can be seen at the end, it still detaches from the logarithmic curve as $x$ increases.

Looking at how the approximation was obtained, you can notice that further differentiation could be used to obtain other approximations for the logarithm function, since for every differentiation I will get the form: $$\mathbb{P}_1(x) \approx \mathbb{P}_2(x)+\ln(x)\cdot \mathbb{P}_3(x) \Rightarrow \ln(x) \approx \frac{\mathbb{P}_1(x)-\mathbb{P}_2(x)}{\mathbb{P}_3(x)}$$

with $\mathbb{P}_i(x)$ some polynomial-alike functions (with roots and fractions of combination of roots and polynomials).

As example, the next step of differentiation leads to the approximation: $$\ln(x)\approx f_2(x) = \frac12 (x-1)^3 \left( \frac{x}{(x-1)^3}-\frac{1}{x(x-1)^3}+\frac{x^{\ln(2)-2}\left(\ln(2)-1-x^{\ln(2)}\right)}{\left(1+x^{\ln(2)}\right)^2}\right)$$

which can be seen in the following figure 3, sacrificing accuracy at the beginning, looks like it is converging for larger values, so I wonder if this process could be indeed giving approximations for the logarithmic function.

2nd approximation


Main question_______________

Considering that the derivative of the logarithm is $\ln'(x) = \frac1x$, that $\frac{d^n}{dx^n}\left(\frac{x}{x-1}\right) = \frac{n!}{(x-1)(1-x)^n}$ and $\frac{d^n}{dx^n}\left(\frac{1}{x}\right) = \frac{(-1)^{n}\ n!}{x^{n+1}}$, and the product rule for higher derivatives $d^{(n)}(uv) =\sum\limits_{k=0}^n {n \choose k}d^{(n-k)}(u)\ d^{(k)}(v)$, if I didn't messed up the nth-order approximation is given by:

$$\ln(x) \approx f_n(x) = \dfrac{ \sum\limits_{k=0}^{n-1}{{n-1}\choose k} \dfrac{d^{k}}{dx^{k}}\left(\dfrac{1}{x}\right)\dfrac{d^{n-1-k}}{dx^{n-1-k}}\left(\dfrac{1}{1+x^{-\ln(2)}}\right) - \sum\limits_{k=1}^{n-1}{{n-1}\choose k} \dfrac{d^k}{dx^k}\left(\dfrac{1}{x}\right)\dfrac{d^{n-1-k}}{dx^{n-1-k}}\left(\dfrac{x}{x-1}\right)}{\dfrac{n!}{(x-1)(1-x)^n}} \label{Eq. 1} \tag{Eq. 1}$$

But here I got stuck:

  1. It is possible to use this for finding a formula for large $x$ values? I was aiming to form a fraction keeping only the higher order polynomials in $n$ and later replace it with $x$, but I couldn't find a formula clearer than this one. Here assuming the formula is right, hope you can review it. I don't know if is possible to do factorization on the $\{n-1-k\}$ terms, as example.
  2. It is converging to the logarithm for large $x$ values as $n$ increase? I don't really know if this is true.

Motivation____________

After receiving a comment criticizing this answer, I realized that approximating $\ln(x) \approx x(x^{\frac1x}-1)$ for large numbers don't going to be useful since looks like $x^{\frac1x}$ is even harder to calculate (even when it looks accurate at large values). So assuming that $x^b$ with $b>0$ is not difficult in this sense (I hope it is not, I am not really sure if make sense the use of $x^{\ln(2)}$), I was aiming to find some ration of polynomials $\ln(x) \approx f_n(x) = \dfrac{\mathbb{P}_1(x^n)}{\mathbb{P}_2(x^n)}$ and later see how behaves replacing $n \leftarrow x$ such $\ln(x) \approx f(x) = \dfrac{\mathbb{P}_1(x^x)}{\mathbb{P}_2(x^x)}$.


Added later

Using Stirling's approximation I could make the following: $$\begin{array}{r c l} \ln(x+1) & = &\ln\left((x+1)!\right)- \ln(x!) \\ \overset{\text{Stirling's}}{\Rightarrow} \ln(x+1)\left(x+\frac12\right) & \approx & 1 + \ln(x)\left(x+\frac12\right) \\ \Rightarrow \ln(x+1)-\ln(x) = \ln\left(1+\frac{1}{x}\right) &\approx & \frac{2}{2x+1} \end{array}$$

But unfortunately I haven't found yet a better simple approximation than $f_1(x)$, neither something useful for large $x$ values.


2nd later comments

The formula \eqref{Eq. 1} was corrected since it was originally wrong. I still without found the formula at $n$, but by mistake I recovered in Wolfra-Alpha other known expansion that works better than Taylor's: $$\ln(x) = \sum\limits_{k=1}^\infty \frac{1}{k}\left(\dfrac{x-1}{x}\right)^k$$ and I am afraid that the second sum of \eqref{Eq. 1} lead to that canceling the logarithm on the approximation "equality", but I still missing with constants matching.


Last update

I give up, after trying much more I have should I found, no after complications because of $\dfrac{d^n}{dx^n} \ln(x)$ broke its form when $n=0$ compared with $n>0$, that the approximation take the form: $$n=0:\qquad \dfrac{\ln(1+x^{\ln(2)})}{\ln(2)}\approx \ln(x)\dfrac{x}{x-1}$$ $$n=1:\qquad \dfrac{1}{x+x^{1-\ln(2)}}\approx \dfrac{1}{x-1}-\dfrac{\ln(x)}{(x-1)^2}\quad \Rightarrow f_1(x)$$

$$\begin{array}{r c l} n>1: \ln(x) & \approx & \sum\limits_{k=1}^n \frac{1}{k}\left(\frac{x-1}{x}\right)^k +\frac{1}{n}\left(\frac{x-1}{x}\right)^n \frac{(x-1)}{(1+x^{\ln(2)})}\sum\limits_{q=0}^{n-1} \sum\limits_{k=0}^q \sum\limits_{j=0}^k \dfrac{(-1)^{j+k-q}{k \choose j}}{(1+x^{-\ln(2)})^k}{(k-j)\ln(2) \choose q} \\ \Rightarrow n\geq 1: \ln(x) & \approx & \sum\limits_{k=1}^{n-1} \frac{1}{k}\left(\frac{x-1}{x}\right)^k + \frac{1}{n}\left(\frac{x-1}{x}\right)^n \left(\dfrac{x+x^{\ln(2)}}{1+x^{\ln(2)}}\right) + \frac{1}{n}\left(\frac{x-1}{x}\right)^n \frac{(x-1)}{(1+x^{\ln(2)})}\sum\limits_{q=1}^{n-1} \sum\limits_{k=0}^q \sum\limits_{j=0}^k \dfrac{(-1)^{j+k-q}{k \choose j}}{(1+x^{-\ln(2)})^k}{(k-j)\ln(2) \choose q} \end{array}$$

Unfortunately I wasn't able to reduce it more neither related it to $f_1(x)$ (which works amazing for small $n$), and also from the graph in Desmos looks like the approximation works good at low $n$ values for $x<1000$ but later it don't gain too much power, but I believe it is converging to the logarithm: the first part do approach the logarithm, and the awful second part, if I am not mistaken, should go to zero for large $n$ and $x$ since are proportional to $\frac{1}{n}$ and because the degree of the polynomials on each term's numerator and denominator matches (hope someone could prove it - I wasn't able to show it).

Another interesting thing is that using the approximation: $$\ln(x) \approx \frac{(x-1)}{x^{\ln(2)}}\ln(1+x^{\ln(2)})$$

and introducing there the classic approximation $\sum_{k=1}^n \frac{1}{k}\left(\frac{x-1}{x}\right)^k$ you get something that arrives a bit faster than the original one: $$\ln(x) \approx \frac{(x-1)}{x\ln(2)}\sum\limits_{k=1}^n \frac{1}{k}\left(\dfrac{x^{\ln(2)}}{1+x^{\ln(2)}}\right)^k$$

5

There are 5 best solutions below

3
On

Answers:

  1. No. Or more precisely, the "drudgery" is not worth it. By drudgery I mean the computational effort. There are simpler procedures.
  2. Maybe yes.

About the computational effort.
IMO, in contrast to the function you suggest the continued fraction Euler shows in §24..25 of Eneström No. 606 is much simpler. I converted the CF to the 6th level, what resulted in $$\log(x)\approx f(x)=\frac{\mathrm{137}x^5+\mathrm{1625}x^4+\mathrm{2000}x^3-\mathrm{2000}x^2-\mathrm{1625}x-\mathrm{137}}{\mathrm{30}\left(x^5+\mathrm{25}x^4+\mathrm{100}x^3+\mathrm{100}x^2+\mathrm{25}x+1\right)}$$ Plotting $f(x)-\log(x)$ unveils disappointing differences for larger and lower $x$ values. f(x)-log(x)
Note: $x=\displaystyle 10^z$
(In addition I found a typo in the original paper and few more in the German translation).

"How is a logarithm computed?" asked Agent Smith in a comment (I guess he knows it) -- if an approximation to 10..12 or few more digits is satisfactory, see HP Journal Vol. 29, No. 8, p. 29 ff.

1
On

First, it should be noted that computers store floating-point numbers in what's essentially scientific notation, with a sign bit ($+$ or $-$), fixed-width, “significand” or “mantissa”, and exponent indicating a power of two:

$$x = (-1)^sm2^p, s \in \{0, 1\}, m \in [1, 2), p \in \mathbb{Z}$$

It's a bit more complicated than that in reality because zero, denormals, infinity, and NaN are special cases. But given that you're taking a logarithm, I'm assuming that you only care about positive finite numbers anyway. For these, we have:

$$\log x = \log m + p\log 2$$

The upshot of this is that if you have a log function that's accurate on the interval $[1, 2]$, you can easily construct one that's accurate $\forall x > 0$. So it's not a problem if your approximation can't handle large $x$'s “directly”.

import math

LOG2 = 0.6931471805599453

def approx_log(x):
    if x <= 0:
        raise ValueError('domain error')
    m, p = math.frexp(x)
    m *= 2; p -= 1 # Adjust for frexp using [1/2, 1) instead of [1, 2)
    log_m = ## some approximation of log(m) that works for 1 <= m < 2
    return log_m + p * LOG2

So let's implement your $f_1$ function.

from math import log

def f1(x):
    y = x ** log(2)
    return (x - 1) * (x + y) / (x * (1 + y))

By evaluating f1 at the values [1 + n/10**6 for n in range(10**6 + 1)], it appears that the maximum error abs(f1(x) - log(x)) is about $0.002075$. Not too bad.

But let's consider alternatives. For example, consider the simple quadratic polynomial:

$$g_2(x) = -0.24021388846266806x^2 + 1.4137888459479493x - 1.1735749574852814$$

(I derived this by interpolating the endpoints $(1, 0)$ and $(2, \ln 2)$, and then experimentally finding a third parameter to minimize the error.) This approximation has an absolute error of approximately $0.005293$. So your $f_1$ is better.

What about a cubic?

$$g_3(x) = 0.1103629380824423x^3 - 0.7345596458883537x^2 + 2.1242855516479104x - 1.500088843841999$$

This one has a maximum error of only $0.000609$, which beats your $f_1$ by a factor of 3.4.

Now, let's consider computation time:

  • $g_3$ can be evaluated using Horner's Method, using 3 multiplications and 3 addition/subtraction operations.
  • $f_1$ as coded above uses 1 exponentiation, 2 multiplications, and 3 addition/subtraction operations.

So, if you subtract what they have in common, it boils down to 1 multiplication version 1 exponentiation.

Now, what are you going to do about that exponentiation? It's not a built-in FPU operator like multiplication is. And you can't rewrite $x^{\ln 2}$ as $e^{\ln 2\ln x}$, because you need an accurate log function in order to do that. Maybe you can write an accurate approximation function for $x^{\ln 2}$, but can it be faster than one multiplication? I doubt it.

In conclusion, your $f_1$ can be beaten, in both accuracy and computation time, by a simple cubic polynomial. So, while your approach may be intriguing, I just don't see it having any real practical use for computing logarithms.

0
On

Hasting's classic "APPROXIMATIONS FOR DIGITAL COMPUTERS" published in 1954 contains a number of nice approximations for $\log_{10}$. Here is the most accurate one.

For $\dfrac1{\sqrt{10}} \le x \le \sqrt{10} $, let $z=\dfrac{x-1}{x+1}$. Then $\log_{10}(x) \approx c_1z+c_3z^3+ c_5z^5+c_7z^7 +c_9z^9$ with

$c_1 = .8685,91718\\ c_3 = .2893,35524\\ c_5 = .1775,22071\\ c_7 = .0943,76476\\ c_9 = .1913,37714 $

with an error less than $1.5*10^{-7}$.

A more recent collection is "Computer Approximations" by Hart, Cheney, Lawson, Maehley, Mesztenyi, Rice, Thacher, and Witzgall, published in 1968. This has a number of quite accurate approximations to log.

0
On

Too long for a comment.

As said in comments, you do not need to cover a wide range of $x$.

What you can do is to use $$\log(x)=-\sum_{n=1}^\infty \frac{(-1)^n}{n}\,(x-1)^n$$ and make it a $[n,n]$ Padé approximant $P_n$ which, after simplifications will write $$P_n=\frac {\sum_{n=0}^\infty a_n \,x^n } {\sum_{n=0}^\infty b_n \,x^n }$$

For example $$P_7=\frac {(x-1) \left(363+10310 x+58673 x^2+101548 x^3+58673 x^4+10310 x^5+363 x^6\right)}{70(1+x)\left(11+48 x+393 x^2+832 x^3+393 x^4+48 x^5+x^6\right) }$$ whose error is $$\frac{(x-1)^{15}}{176679360}$$

For the useful range already mentioned by @marty cohen, the maximum error is $3.10\times 10^{-8}$.

7
On

From what I got you would like to find an approximation of $\ln(x)$ for big values of $x$.

If so then you can consider $$\ln(x)\approx\frac{\pi}{2M\big(1,2^{2-m}/x\big)}-m\ln(2)$$ where $M(x,y)$ is arithmetic-geometric mean.

Here are $\log$ for first $10$ natural numbers calculated directly using Wolfram Mathematica:

Table[N[ Log[x], 20], {x, 2, 10}]

$${0.69314718055994530942, 1.0986122886681096914, \ 1.3862943611198906188, 1.6094379124341003746, 1.7917594692280550008, \ 1.9459101490553133051, 2.0794415416798359283, 2.1972245773362193828, \ 2.3025850929940456840}$$

and using arithmetic-geometric mean:

Table[N[Pi/(2 ArithmeticGeometricMean[1, 2^(2 - m)/x]) - m Log[2], 
   20], {x, 2, 10}] //. m -> 10

$${0.693147180559945, 1.098612288668110, 1.386294361119891, \ 1.609437912434100, 1.791759469228055, 1.945910149055313, \ 2.079441541679836, 2.197224577336219, 2.302585092994046}$$