Post Reply 
Estimate logarithm quickly
08-21-2021, 03:39 PM
Post: #1
Estimate logarithm quickly
Estimate logorithm quickly, with tiny errors.
source: Dead Reckoning: Calculating without Instruments, p. 123

n*n = (n+x) * (n-x) * (n*n/(n*n-x*x))

2*ln(n) = ln(n+x) + ln(n-x) + ln(n*n/(n*n-x*x))

ln(n*n/(n*n-x*x)) = ln(1-x*x/(n*n-x*x)) = 2*atanh(x*x/(2*n*n-x*x))

Let a = ln(n+x), b = ln(n-x), c = n/x, and halved both side.

ln(n) = (a+b)/2 + atanh(1/(2*c*c-1))

atanh(1/(2*c*c-1)) > atanh(1/c) / (2c) = (a-b)/2 / (2c)

ln(n) > (a+b)/2 + (a-b)/2 / (2c)

If we divide both side by ln(base), inequality still holds
We can thus generalize for log of any base.

Example, n=3, x=1

a = log(3+1) = 2*log(2)
b = log(3-1) = log(2)
c = 3/1 = 3

log(3) > 3*log(2)/2 + (log(2)/2) / (2*3) = 19/12 * log(2)
log(3)/log(2) > 19/12

Example, n=3, x=1.5

a = log(3+1.5) = 2*log(3)-log(2)
b = log(3-1.5) = log(3)-log(2)
c = 3/1.5 = 2

log(3) > (3*log(3)-2*log(2))/2 + (log(3))/2 / (2*2) = -log(2) + 13/8*log(3)
log(2) > 5/8*log(3)
log(3)/log(2) < 8/5

19/12 < log2(3) < 8/5

Both extreme happened to be convergents of log2(3)

C:\> spigot -C log(3)/log(2) | head
1/1
2/1
3/2
8/5
19/12
65/41
84/53
485/306
1054/665
24727/15601
Find all posts by this user
Quote this message in a reply
08-21-2021, 03:58 PM
Post: #2
RE: Estimate logarithm quickly
(08-21-2021 03:39 PM)Albert Chan Wrote:  atanh(1/(2*c*c-1)) > atanh(1/c) / (2c) = (a-b)/2 / (2c)

This is not yet proofed (n-x > 0, c = n/x > 1)
To proof it, it required this inequality (also not proofed, suggestions welcome):

(1-w)^(2-w) * (1+w)^(2+w) < 1,       where 0 < w < 1

XCAS> series((1-w)^(2-w) * (1+w)^(2+w), w, 0, 10)

1 - w^4/3 - 4*w^6/15 - 10*w^8/63 - 4*w^10/45 + w^11*order_size(w)

Taylor series suggested the expression peaked at w = 0, and plots support it.
For now, assuming it is true ...

(1-w)^(2-w) * (1+w)^(2+w) < 1
(1-w*w)^2 < ((1-w)/(1+w))^w
(1/(1-w*w))^2 > ((1+w)/(1-w))^w

Let c = 1/w, c > 1

(c*c/(c*c-1))^2 > ((1+1/c)/(1-1/c)) ^ (1/c)

ln(c*c/(c*c-1)) * 2 > ln((1+1/c)/(1-1/c)) / c

Simpilfy with identity ln(z) = 2*atanh((z-1)/(z+1))

(2*atanh(1/(2*c*c-1)) * 2 > (2*atanh(1/c)) / c

atanh(1/(2*c*c-1) > atanh(1/c) / (2c)
Find all posts by this user
Quote this message in a reply
08-21-2021, 05:41 PM
Post: #3
RE: Estimate logarithm quickly
(08-21-2021 03:58 PM)Albert Chan Wrote:  To proof it, it required this inequality (also not proofed, suggestions welcome):

(1-w)^(2-w) * (1+w)^(2+w) < 1,       where 0 < w < 1

My attempt to proof this, by showing there is no local extremum for 0 < w < 1
To get the extrema of above is messy. Instead, we solve extrema of ln(expression)

XCAS> f := (2-w)*ln(1-w) + (2+w)*ln(1+w)
XCAS> diff(f,w)

-ln(1-w) - (2-w)/(1-w) + ln(1+w) + (2+w)/(1+w)

XCAS> numer(ans()) // local extrema if df/dw = 0

-ln(w+1) + ln(-w+1) + w^2*ln(w+1) - w^2*ln(-w+1) + 2*w

Divide all by 2, then replace with identity atanh(x) = ln((x+1)/(1-x))/2

(w^2-1) * atanh(w) + w = (w^2-1) * (w + w^3/3 + w^5/5 + ...) + w

All coefficients of w's is positive, 1/(2k-1) - 1/(2k+1) = 2/(4k^2-1) > 0

XCAS> series((w^2-1)*atanh(w) + w, w,0,10)

2*w^3/3 + 2*w^5/15 + 2*w^7/35 + 2*w^9/63 + w^11*order_size(w)

Within 0 < w < 1, there is no local extremum, f is bounded by edges.

XCAS> limit(f, w=0, +1)      → 0
XCAS> limif(f, w=1, -1)       → -infinity

f = ln((1-w)^(2-w) * (1+w)^(2+w)) < 0       → (1-w)^(2-w) * (1+w)^(2+w) < 1
Find all posts by this user
Quote this message in a reply
08-21-2021, 06:31 PM
Post: #4
RE: Estimate logarithm quickly
Albert, thanks for the book reference.
(Surprisingly, that's one I don't have...yet.)

ENTER > =
Find all posts by this user
Quote this message in a reply
08-21-2021, 11:42 PM
Post: #5
RE: Estimate logarithm quickly
This is an amazing formula, for estimating ln(n) to 5 sig. digits, n = [1/√2, √2]
(again, from "Dead Reckoning", p. 132)

ln(n) ≈ 6*(n-1) / (4*√n + n+1)

use 1/ln(10) = 0.4342945 ... to scale back for log10.

Example, for log10(3):

3^2 = 0.9e1

log10(0.9) ≈ 0.4343 * 6 * -.1 / (4*√.9 + 1.9) ≈ -0.04576
log10(3) ≈ (-0.04576 + 1) / 2 = 0.47712


Or, scale by 2^k factor, log10(2) ≈ 0.30103

3 * 4 = 1.2e1

log10(1.2) = 0.4343 * 6*.2/(4*√1.2 + 2.2) ≈ 0.07918
log10(3) ≈ 0.07918 + 1 - 2*0.30103 = 0.47712
Find all posts by this user
Quote this message in a reply
08-23-2021, 06:44 AM
Post: #6
RE: Estimate logarithm quickly
(08-21-2021 03:39 PM)Albert Chan Wrote:  Estimate logorithm quickly, with tiny errors.
...
Both extreme happened to be convergents of log2(3)

C:\> spigot -C log(3)/log(2) | head
1/1
2/1
3/2
8/5
19/12
65/41
...

Just to note, this is Simon Tatham's arbitrary-precision spigot calculator, which deals in continued fractions as well as decimals:
https://www.chiark.greenend.org.uk/~sgta...html#cfrac
Find all posts by this user
Quote this message in a reply
10-06-2021, 10:58 PM
Post: #7
RE: Estimate logarithm quickly
(10-06-2021 05:12 PM)Eric Rechlin Wrote:  Computing Logarithms by Integration - Richard Schwartz (17:02)

Thanks for the video. I finally "get" Doerfler's formula Smile

ln(n) = (n-1) / ∫(n^t, t, 0, 1) ≈ (n-1) / ((1 + 4*√n + n)/6) // Simpson's Rule.

I also learned about Boole's Rule.
It is really just Richardson extrapolation, from Simpson's Rule.

CAS> S1 := [1,0,4,0,1] ./ 6
CAS> S2 := [1,4,2,4,1] ./ 12
CAS> B1 := S2 + (S2-S1)/(4^2-1)     → [7/90, 16/45, 2/15, 16/45, 7/90]
CAS> horner(B1, e^0.25) / (e-1)       → 1.00000050019 // estimated 1/ln(e), Boole's Rule

This should match extrapolation of 2 Simpson's.
Note: we work with reciprocal log, to simplify Richardson's extrapolation.

CAS> D(x) := horner([1,4,1], x^0.5)/(6*(x-1)) // Doerfler's formula, for 1/ln(x)
CAS> s1 := D(e)             → 1.00033715273
CAS> s2 := D(e^0.5)/2   → 1.00002154097
CAS> horner([1,4,2,4,1], e^0.25)/(12*(e-1))   → 1.00002154097

Note: s2 is equivalent to more messy [1,4,2,4,1] Simpson's. version.
Thus, it is a good idea to reduce log arguments, close to 1

Extrapolate for Boole's Rule result:

CAS> s2 + (s2-s1)/(4^2-1)              → 1.00000050019
Find all posts by this user
Quote this message in a reply
10-07-2021, 06:25 PM
Post: #8
RE: Estimate logarithm quickly
I have a very easy to remember method if I want to sketch log/log or log/lin data with hand:
log2 --> 0.3
log3 --> 0.5
log5 --> 0.7

These values approx equally spaced on the log scale and with log1 and log10 values there are 5 points on one decade.

Csaba
Find all posts by this user
Quote this message in a reply
10-18-2021, 01:02 PM (This post was last modified: 10-18-2021 04:30 PM by Albert Chan.)
Post: #9
RE: Estimate logarithm quickly
(10-06-2021 10:58 PM)Albert Chan Wrote:  
(10-06-2021 05:12 PM)Eric Rechlin Wrote:  Computing Logarithms by Integration - Richard Schwartz (17:02)

Thanks for the video. I finally "get" Doerfler's formula Smile

Both methods actually are one and the same.
Doerfler's Borchardt's algorithm + Richardson extrapolation: d(k,k) ≈ ∫(n^t,t,0,1) = (n-1)/ln(n)

CAS> d00 := a0 := (1+n)/2

CAS> g0 := sqrt(n)
CAS> d01 := factor(a1 := (a0+g0)/2)       → (n + 2*n^(1/2) + 1)/4
CAS> d11 := factor(d01 + (d01-d00)/3)    → (n + 4*n^(1/2) + 1)/6

CAS> g1 := sqrt(a1*g0)       → (sqrt(sqrt(n)*(n+1)+2*n))/2

This seemingly looking "AGM" is really to evaluate doubled points. (*)
We can rewrite g1 = (√√n) * (1+√n)/2, expanded:

CAS> g1 := (n^(1/4) + n^(3/4))/2
CAS> d02 := factor(a2 := (a1+g1)/2)       → (n + 2*n^(3/4) + 2*n^(1/2) + 2*n^(1/4) + 1)/8
CAS> d12 := factor(d02 + (d02-d01)/3)    → (n + 4*n^(3/4) + 2*n^(1/2) + 4*n^(1/4) + 1)/12
CAS> d22 := factor(d12 + (d12-d11)/15)   → (7*n + 32*n^(3/4) + 12*n^(1/2) + 32*n^(1/4) + 7)/90

(*) We can "see" doubling of points by thinking weights as a "number"
weight of a3 = 122222221 = 11111111*11
weight of g2 = 1010101 = 11111111/11
weight of g3 = √(a3*g2) = 11111111
Find all posts by this user
Quote this message in a reply
Post Reply 




User(s) browsing this thread: 2 Guest(s)