Estimate logarithm quickly

[Albert Chan]

Instead of Romberg integration, we can also do Gaussian Quadrature, to estimate ln(n)

Shifting points and weights to limit of 0 to 1, for 3 points, we have:

Code:

point         weight
0.5            8/18
0.5 ± √(0.15)  5/18

ln(n) ≈ (n-1) * 18/(((n^p + n^-p)*5 + 8)*√n)       , where p = √(0.15)

Above formula give ln(e) ≈ 1.0000004796, slightly better than Booles' rule.
However, it maybe expensive to calculate n^p, n^-p

The good news is (n^p + n^-p) taylor series, half the terms get cancelled.
Only even powers of p remains, no square roots anywhere.

CAS> s := (1+x)^p + (1+x)^-p
CAS> series(s)

2 + p^2*x^2 - p^2*x^3 + (1/12*p^4+11/12*p^2)*x^4 + (-1/6*p^4-5/6*p^2)*x^5 + x^6*order_size(x)

CAS> pade(s, x, 4, 3) (p=sqrt(0.15))

(-3.5*x^2-24*x-24)/(-0.85*x^2-12*x-12)

Replacing s with above pade approximation, we have:

ln(n) ≈ g - g/(2.7 + 24/g^2)       , where g = (n-1)/√n (*)

With this, we have log(e) ≈ 1.00016029872
Not as good, but still twice as good as Doerfler's formula (log(e) ≈ 1/1.00033715273)

If we limit n = [√0.5, √2], this Pade formula is at least 17 times better.

lua> n = 1.2586
lua> (n-1) * 6/(1 + 4*sqrt(n) + n) -- Doerfler formula
0.22999976897818544
lua> g = (n-1)/sqrt(n)
lua> g - g / (2.7 + 24/g^2) -- Pade formula
0.22999999788822215
lua> log(n)
0.2299999921106961


(*) If we consider whole expression (not just 3 points), we still get the same result.

CAS> pade(ln(1+x)*sqrt(1+x)/x,x,4,3)       → (17*x^2+240*x+240)/(27*x^2+240*x+240)

→ ln(1+x)
≈ x/sqrt(1+x) * (1 - 1/(27/10 + 24*(1+x)/x^2))
≈ x - x^2/2 + x^3/3 - x^4/4 + x^5/5 - x^6/6 + x^7/6.992 - x^8/7.962 + x^9/8.897 - x^10/9.786 + ...