2

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$$

Joako
  • 1,380
  • 1
    Neat! Question: How is a logarithm computed? – Agent Smith Mar 20 '24 at 23:44
  • @AgentSmith I don't really know, I believe that through Taylor's ¿?, I am using Wolfram-Alpha and if I am not mistaken they use Meiger's G functions (which I don't understand). I am assuming that a $x^{\ln(2)}$ root is not something difficult to calculate, but I am not sure neither (I think is just a polynomial, hope I am not mistaken). I start by exploring the formula I found, and get the doubt because someone make me note $\ln(x) \approx x x^{1/x}-x$ is not make too much sense, which I found by changing $2^n$ to $x$ from this another answer. – Joako Mar 20 '24 at 23:51
  • Well, my hunch is, from a few approximations I've seen (e.g. Brahmagupta's formula for sine), the error in the approximating formula is least under $2$ conditions: a) when $x$ is really, really small and b) when $x$ is really, really big. Most of the deviations seem to occur at "medium" $x$ values. Is there a way to fix that? Je ne sais pas. – Agent Smith Mar 20 '24 at 23:58
  • 4
    Your approach is really interesting, but it doesn't necessarily need to extend to large values by itself. Once you have a way to compute $\ln$ for small values, you can use the fact that $\ln(a^n b)=n \ln a +\ln b$ to take a large value, scale it down by $a^n$, and then compute the logarithm of that smaller scaled value. – Alex K Mar 21 '24 at 00:58
  • @AlexK that is interesting, I didn't think on that... I will try to see if I can use it – Joako Mar 21 '24 at 02:38
  • 1
    @Joako R. P. Brent's MP from the 1970s computes log1p(x) for x small in magnitude with Newton iteration based on expm1(). expm1() itself can be every efficiently computed using the Taylor expansion of sinh(), but I do not recall whether this is how MP computes it, or whether it uses a different method. MP computes log(x) for larger x by first computing a rough approximation y=log(x). It then computes exp(y) accurately and finally computes log(x/exp(y)) utilizing the log1p() computation mentioned earlier. – njuffa Mar 21 '24 at 05:06
  • 1
    Euler presents in a paper with Eneström No 606 the approximation of $\log$ by continued fraction. See §24 in a. m. reference. – m-stgt Mar 21 '24 at 22:00
  • @AlexK Unfortunately I have not been able of finding something useful. But I realize, unrelated, that $f(x)+f(x^{-1})=0,\ f(1)=0$ is fulfilled by the logarithm, but also by the approximation $f_1(x)$. This means are many solutions to this functional equation, maybe it could be useful for finding other approximations. – Joako Mar 21 '24 at 22:33
  • 1
    Computer floating point numbers are effectively represented in scientific notation: $x = m2^p$, where $m \in [1, 2)$. You can thus compute $\log x = \log m + p\log 2$, and thus having a log function that's accurate on the interval $[1, 2]$ automatically gives you a log function for all finite positive numbers. Thus, there's no need to "directly" computer the log of a large number. – Dan Mar 22 '24 at 05:35
  • @Dan How you could use it with the $f_1(x)$ example? it shows to be quite accurate in your mentioned range $1,\ 2]$, but it get really bad as $x$ go big ($x>50$ as example) – Joako Mar 22 '24 at 15:45
  • @Joako: How exactly are you intending to compute $x^{\ln(2)}$? The conventional way is $e^{\ln(2)\ln(x)}$, but you need a log function before you can do that. – Dan Mar 22 '24 at 16:11
  • @Dan I am assuming is just a root, like taking $x^{3/4}$ as example, but I stated at the end of my question, I am not sure about if this is true. – Joako Mar 22 '24 at 17:57
  • @TurlocTheRed It looks like works only for $x\in [1,\ 3/2]$ , Maybe there is a typo? – Joako Mar 22 '24 at 22:15

5 Answers5

3

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.

Dan
  • 14,978
  • Thanks for taking the time of giving a detailed answer. I understand your point about the calculation procedure now, but when I am asking about "big numbers" I am using it in the sense of the Stirling's approximation, and I am not sure if it is in the same line of your answer. By the way, since you focused your question in the floating point description, a few months ago I accidentally found a function where floating point description get screwed, at least on graphing tools. – Joako Mar 22 '24 at 21:22
2

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.

m-stgt
  • 301
  • Thanks for taking the time to answer. The link to the HP Journal was really interesting. – Joako Mar 22 '24 at 15:42
  • I have compared here the answer you shared which have polynomials of order $5$, and it looks it works worst that the function $f_1(x)$ which is the simplest example I show which have order $< 2$, so I don't think my approach is more complicated, since looks like works better using lower order polynomials. – Joako Mar 22 '24 at 15:42
  • 1
    @Joako -- Yes, you are correct, I was astonished that Euler's CF is useful in a quite limited input range only. That is why I state "disappointing differences for larger and lower $x$ values." But, as Dan outlined in his comment you may "normalize" input to this acceptable range of accuracy, $\displaystyle \log m2^p = \log m + p\log 2$. Analogical did once HP in their pocket calculators based on BCD. – m-stgt Mar 22 '24 at 16:19
1

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}$.

1

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}$$

  • Thanks for answering. It is an interesting result, unfortunately I couldn't test it for really large values in Wolfram-Alpha since looks like the built-in functions got out of range. – Joako Mar 29 '24 at 04:41
  • By the way, it is not answering what I am asking for. I think it could be more suitable for this question. – Joako Mar 29 '24 at 04:43
  • How you get this approximation? unfortunately I cannot test it for $x>1500$, but below it looks converges faster to $\ln(x)$ than $x(x^{1/x}-1)$ – Joako Mar 29 '24 at 20:00
  • For example for $x=2000$ and $m=20$ the calculated result is $7.60090245954208238$ and the exact logarithm is $7.6009024595420823615$ – Gevorg Hmayakyan Mar 29 '24 at 21:05
  • $M$ is easy to calculate, there is a few methods how to efficiently calculate it. Especially there is an iterative method and elliptic integral method. – Gevorg Hmayakyan Mar 29 '24 at 21:07
  • I mean how you got this formula? I guess that there is some background calculations that explains its success on fitting the logarithm – Joako Mar 29 '24 at 21:18
  • 1
    Please look at this: https://en.wikipedia.org/wiki/Logarithm#Arithmetic%E2%80%93geometric_mean_approximation – Gevorg Hmayakyan Mar 29 '24 at 21:23
0

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.

marty cohen
  • 107,799