Natural logarithm of 2 [HP-15C/HP-42S/Free42 & others]

[Albert Chan]

(11-11-2019 06:14 PM)Gerson W. Barbosa Wrote:  n = 1663 will be enough for 1000 digits:

It may be better to make n even, avoiding 1 square root:

\(\Large \frac{\pi}{\log 2}\normalsize ≈ AGM(2n, \frac{2n}{2^{n-2}}) = AGM(n + \frac{n}{2^{n-2}}, \frac{2n}{2^{n/2-1}}) \)

With n=1662, loop count down to 17.

Code:

OPTION ARITHMETIC DECIMAL_HIGH
LET n = 1662
LET nd = 1000
PRINT "Log(2) = "
LET t = TIME 
LET c = 2^(n/2-1)
LET a = n + n/(c*c)
LET b = (n+n)/c
FOR i = 1 TO 17
   LET c = a
   LET a = (a + b)/2
   LET b = SQR(b*c)
NEXT i
LET c = (a + b)/2
LET a = c - (a-b)^2/(16*c)
LET g = PI/a
LET r = TIME - t
LET r$ = STR$(g)                           
PRINT
PRINT "0";
PRINT r$(0:1);
FOR i = 2 TO nd + 1
   PRINT r$(i:i);
   IF MOD((i - 1),10) = 0 THEN PRINT " ";
   IF MOD((i - 1),50) = 0 THEN 
      PRINT
      PRINT "  ";
   END IF
NEXT i
IF MOD (i - 2,50) <> 0  OR nd = 0 THEN PRINT
PRINT 
PRINT "Runtime: ";
PRINT  USING "0.##": r;
PRINT " seconds"
END

[Paul Dale]

(11-13-2019 05:40 PM)Gerson W. Barbosa Wrote:  Former forum member Mike (Stgt) sent me two RPN programs which give perfect 10-digit and 12-digit results instantaneously on RPN calculators...

Thanks for this, I'd not found enough time to have a proper go.


Pauli

[Gerson W. Barbosa]

(11-14-2019 12:04 AM)Albert Chan Wrote:  

(11-11-2019 06:14 PM)Gerson W. Barbosa Wrote:  n = 1663 will be enough for 1000 digits:

It may be better to make n even, avoiding 1 square root:

\(\Large \frac{\pi}{\log 2}\normalsize ≈ AGM(2n, \frac{2n}{2^{n-2}}) = AGM(n + \frac{n}{2^{n-2}}, \frac{2n}{2^{n/2-1}}) \)

With n=1662, loop count down to 17.

Code:

OPTION ARITHMETIC DECIMAL_HIGH
LET n = 1662
LET nd = 1000
PRINT "Log(2) = "
LET t = TIME 
LET c = 2^(n/2-1)
LET a = n + n/(c*c)
LET b = (n+n)/c
FOR i = 1 TO 17
   LET c = a
   LET a = (a + b)/2
   LET b = SQR(b*c)
NEXT i
LET c = (a + b)/2
LET a = c - (a-b)^2/(16*c)
LET g = PI/a
LET r = TIME - t
LET r$ = STR$(g)                           
PRINT
PRINT "0";
PRINT r$(0:1);
FOR i = 2 TO nd + 1
   PRINT r$(i:i);
   IF MOD((i - 1),10) = 0 THEN PRINT " ";
   IF MOD((i - 1),50) = 0 THEN 
      PRINT
      PRINT "  ";
   END IF
NEXT i
IF MOD (i - 2,50) <> 0  OR nd = 0 THEN PRINT
PRINT 
PRINT "Runtime: ";
PRINT  USING "0.##": r;
PRINT " seconds"
END

I`ve just tested your optmization. For conciseness, here are the last lines of the ouput:

...
7891923723 1467232172 0534016492 5687274778 2344535347

Runtime: 0.03 seconds

The running times vary, depending on whatever other tasks the OS is doing. The ones I've posted are the minimum values I've observed for each example.

This appears to be the record here for one thousand correct digits.

[Albert Chan]

Another way to speedup log(2) calculations is to use ... uh ... log(2)

Let r = roughly the value of log(2)

log2 = r + log2 - r = r + log(2 exp(-r))

Last term should be tiny:

log(2 exp(-r)) = log(1 + ε) = 2 atanh(x)

With r = 0.6931471805599453, x = ε/(ε+2) ≈ 4.7e-18

→ atanh(x) converge very quickly, +35 decimal digits per loop.

Code:

OPTION ARITHMETIC DECIMAL_HIGH
LET nd = 1000
PRINT "Log(2) ="
LET t = TIME
LET r = 0.6931471805599453    !'log2 16 digits
LET a = r*r
LET x = 1
LET y = 0.5 - r               !'y = exp(-r) - 1/2
FOR k = 3 TO 425 STEP 2       !'k loops = 212
   LET x = x * a / (k*k-k)
   LET y = y + x * (k-r)      !'2 terms at a time
NEXT k
LET x = y / (y+1)             !'log(1+2y) = 2*atanh(x)
LET a = x*x
LET s = x
FOR d = 3 TO 57 STEP 2        !'d loops = 28
   LET x = x*a
   LET s = s + x / d
NEXT d
LET s = r + 2*s               !log2 = r + 2*atanh(x)
LET r = TIME - t

LET r$ = STR$(s)                           
PRINT
PRINT "0";
PRINT r$(0:1);
FOR i = 2 TO nd + 1
   PRINT r$(i:i);
   IF MOD((i - 1),10) = 0 THEN PRINT " ";
   IF MOD((i - 1),50) = 0 THEN 
      PRINT
      PRINT "  ";
   END IF
NEXT i
IF MOD (i - 2,50) <> 0  OR nd = 0 THEN PRINT
PRINT 
PRINT "Runtime: ";
PRINT  USING "0.##": r;
PRINT " seconds"
END

[Gerson W. Barbosa]

(11-14-2019 05:34 PM)Albert Chan Wrote:  

Code:

OPTION ARITHMETIC DECIMAL_HIGH
LET nd = 1000
PRINT "Log(2) ="
LET t = TIME
LET r = 0.6931471805599453    !'log2 16 digits
LET a = r*r
LET x = 1
LET y = 0.5 - r               !'y = exp(-r) - 1/2
FOR k = 3 TO 425 STEP 2       !'k loops = 212
   LET x = x * a / (k*k-k)
   LET y = y + x * (k-r)      !'2 terms at a time
NEXT k
LET x = y / (y+1)             !'log(1+2y) = 2*atanh(x)
LET a = x*x
LET s = x
FOR d = 3 TO 57 STEP 2        !'d loops = 28
   LET x = x*a
   LET s = s + x / d
NEXT d
LET s = r + 2*s               !log2 = r + 2*atanh(x)
LET r = TIME - t

LET r$ = STR$(s)                           
PRINT
PRINT "0";
PRINT r$(0:1);
FOR i = 2 TO nd + 1
   PRINT r$(i:i);
   IF MOD((i - 1),10) = 0 THEN PRINT " ";
   IF MOD((i - 1),50) = 0 THEN 
      PRINT
      PRINT "  ";
   END IF
NEXT i
IF MOD (i - 2,50) <> 0  OR nd = 0 THEN PRINT
PRINT 
PRINT "Runtime: ";
PRINT  USING "0.##": r;
PRINT " seconds"
END

In decimal BASIC you can also use

LET r = LOG(2)               !'log2 16 digits

Either way, the output is

Log(2) =

0.6931471805 5994530941 7232121458 1765680755 0013436025 
  5254120680 0094933936 2196969471 5605863326 9964186875 
  4200148102 0570685733 6855202357 5813055703 2670751635 
  0759619307 2757082837 1435190307 0386238916 7347112335 
  0115364497 9552391204 7517268157 4932065155 5247341395 
  2588295045 3007095326 3666426541 0423915781 4952043740 
  4303855008 0194417064 1671518644 7128399681 7178454695 
  7026271631 0645461502 5720740248 1637773389 6385506952 
  6066834113 7273873722 9289564935 4702576265 2098859693 
  2019650585 5476470330 6793654432 5476327449 5125040606 
  9438147104 6899465062 2016772042 4524529612 6879465461 
  9316517468 1392672504 1038025462 5965686914 4192871608 
  2938031727 1436778265 4877566485 0856740776 4845146443 
  9940461422 6031930967 3540257444 6070308096 0850474866 
  3852313818 1676751438 6674766478 9088143714 1985494231 
  5199735488 0375165861 2753529166 1000710535 5824987941 
  4729509293 1138971559 9820565439 2871700072 1808576102 
  5236889213 2449713893 2037843935 3088774825 9701715591 
  0708823683 6275898425 8918535302 4363421436 7061189236 
  7891923723 1467232172 0534016492 5687274778 2344535347 
  
Runtime: 0.02 seconds

[toml_12953]

(11-14-2019 05:34 PM)Albert Chan Wrote:  Another way to speedup log(2) calculations is to use ... uh ... log(2)

Let r = roughly the value of log(2)

It's not good practice to reuse a variable other than loop indices. Since r is already used for the run time, you should use some variable other than r for LOG(2).

[Albert Chan]

(11-15-2019 10:58 AM)toml_12953 Wrote:  It's not good practice to reuse a variable other than loop indices.
Since r is already used for the run time, you should use some variable other than r for LOG(2).

This faster version renamed variable r as g.

To speedup exp(-g): exp(-g) = exp(-g/n) ^ n, where n = 2^24

This reduce exp() loop count from 212 down to 55.

Code:

OPTION ARITHMETIC DECIMAL_HIGH
LET n = 2^24
LET nd = 1000
PRINT "Log(2) ="
LET t = TIME
LET g = 0.6931471805599453    !'log2 16 digits
LET h = g / n
LET a = h*h
LET x = 1
LET s = 1 - h                 !'s = exp(-h)
FOR k = 3 TO 111 STEP 2       !'loops = 55
   LET x = x * a / (k*k-k)
   LET s = s + x * (k-h)      !'2 terms at a time
NEXT k
LET s = s^n - 0.5             !'s = exp(-g) - 0.5
LET x = s / (s+1)             !'log(1+2s) = 2*atanh(x)
LET a = x*x
LET s = x
FOR k = 3 TO 57 STEP 2        !'loops = 28
   LET x = x*a
   LET s = s + x / k
NEXT k
LET s = g + 2*s               !log2 = g + 2*atanh(x)

LET r = TIME - t
LET r$ = STR$(s)                           
PRINT
PRINT "0";
PRINT r$(0:1);
FOR i = 2 TO nd + 1
   PRINT r$(i:i);
   IF MOD((i - 1),10) = 0 THEN PRINT " ";
   IF MOD((i - 1),50) = 0 THEN 
      PRINT
      PRINT "  ";
   END IF
NEXT i
IF MOD (i - 2,50) <> 0  OR nd = 0 THEN PRINT
PRINT 
PRINT "Runtime: ";
PRINT  USING "0.##": r;
PRINT " seconds"
END

[Gerson W. Barbosa]

(11-17-2019 03:01 PM)Albert Chan Wrote:  To speedup exp(-g): exp(-g) = exp(-g/n) ^ n, where n = 2^24

This reduce exp() loop count from 212 down to 55.

Code:

OPTION ARITHMETIC DECIMAL_HIGH
LET n = 2^24
LET nd = 1000
PRINT "Log(2) ="
LET t = TIME
LET g = 0.6931471805599453    !'log2 16 digits
LET h = g / n
LET a = h*h
LET x = 1
LET s = 1 - h                 !'s = exp(-h)
FOR k = 3 TO 111 STEP 2       !'loops = 55
   LET x = x * a / (k*k-k)
   LET s = s + x * (k-h)      !'2 terms at a time
NEXT k
LET s = s^n - 0.5             !'s = exp(-g) - 0.5
LET x = s / (s+1)             !'log(1+2s) = 2*atanh(x)
LET a = x*x
LET s = x
FOR k = 3 TO 57 STEP 2        !'loops = 28
   LET x = x*a
   LET s = s + x / k
NEXT k
LET s = g + 2*s               !log2 = g + 2*atanh(x)

LET r = TIME - t
LET r$ = STR$(s)                           
PRINT
PRINT "0";
PRINT r$(0:1);
FOR i = 2 TO nd + 1
   PRINT r$(i:i);
   IF MOD((i - 1),10) = 0 THEN PRINT " ";
   IF MOD((i - 1),50) = 0 THEN 
      PRINT
      PRINT "  ";
   END IF
NEXT i
IF MOD (i - 2,50) <> 0  OR nd = 0 THEN PRINT
PRINT 
PRINT "Runtime: ";
PRINT  USING "0.##": r;
PRINT " seconds"
END

Output:

Log(2) =

0.6931471805 5994530941 7232121458 1765680755 0013436025 
  5254120680 0094933936 2196969471 5605863326 9964186875 
  4200148102 0570685733 6855202357 5813055703 2670751635 
  0759619307 2757082837 1435190307 0386238916 7347112335 
  0115364497 9552391204 7517268157 4932065155 5247341395 
  2588295045 3007095326 3666426541 0423915781 4952043740 
  4303855008 0194417064 1671518644 7128399681 7178454695 
  7026271631 0645461502 5720740248 1637773389 6385506952 
  6066834113 7273873722 9289564935 4702576265 2098859693 
  2019650585 5476470330 6793654432 5476327449 5125040606 
  9438147104 6899465062 2016772042 4524529612 6879465461 
  9316517468 1392672504 1038025462 5965686914 4192871608 
  2938031727 1436778265 4877566485 0856740776 4845146443 
  9940461422 6031930967 3540257444 6070308096 0850474866 
  3852313818 1676751438 6674766478 9088143714 1985494231 
  5199735488 0375165861 2753529166 1000710535 5824987941 
  4729509293 1138971559 9820565439 2871700072 1808576102 
  5236889213 2449713893 2037843935 3088774825 9701715591 
  0708823683 6275898425 8918535302 4363421436 7061189236 
  7891923723 1467232172 0534016492 5687274778 2344535347 
  
Runtime: 0.00 seconds

[toml_12953]

(11-17-2019 03:01 PM)Albert Chan Wrote:  This faster version renamed variable r as g.

To speedup exp(-g): exp(-g) = exp(-g/n) ^ n, where n = 2^24

This reduce exp() loop count from 212 down to 55.

Wow! That's some reduction. Good job!

[Albert Chan]

Inspired by Valentin Albillo post on estimating Pi using random numbers.
And, the proof for it: https://math.stackexchange.com/questions...379#885379

LN(2) = 2 * probability of integer part of RND/RND is odd

10 INPUT K @ N=0 @ FOR I=1 TO K @ N=N+MOD(IP(RND/RND),2) @ NEXT I @ DISP 2*N/K @ GOTO 10

> RANDOMIZE 1
> RUN
? 1E3     → .68
? 1E3     → .756
? 1E4     → .6806
? 1E4     → .6938

[EdS2]

(03-23-2020 04:34 PM)Albert Chan Wrote:  Inspired by Valentin Albillo post on estimating Pi using random numbers.
And, the proof for it: https://math.stackexchange.com/questions...379#885379

Oh, very many thanks for that link! A simple proof and exactly the graphic that I've been wanting to see!

Quote:LN(2) = 2 * probability of integer part of RND/RND is odd

Goodness me!

[Albert Chan]

(03-23-2020 04:34 PM)Albert Chan Wrote:  LN(2) = 2 * probability of integer part of RND/RND is odd

We can improve LN(2) estimate by scaling RND/RND
Let p(m) = probability of integer part of m*RND/RND is odd

XCas> [log2] := expand(solve(p(m) = floor(m/2)/(2*m) + (m/2)*(log2 + c), log2))     → [2*p(m)/m-floor(m/2)/m^2-c]
XCas> c := (m>1)*sum((-1)^k/k, k=1..2*floor(m/2))

XCas> expand(subst(log2, m=1))     → 2*p(1)
XCas> expand(subst(log2, m=2))     → p(2)+1/4
XCas> expand(subst(log2, m=3))     → 2*p(3)/3+7/18
XCas> expand(subst(log2, m=10))   → p(10)/5+1501/2520

> 10 DEF FNL(K) @ N=0
> 20 FOR L1=1 TO K @ N=N+MOD(IP(10*RND/RND),2) @ NEXT L1
> 30 FNL=N/(5*K)+1501/2520 @ END DEF
> RUN

> FIX 5 @ RANDOMIZE 1
> FOR I=1 to 5 @ FNL(1e4) @ NEXT I
0.69063
0.69441
0.69329
0.69347
0.69377

[Thomas Klemm]

(10-12-2019 12:32 AM)Gerson W. Barbosa Wrote:  By fast I mean under 10 seconds on the HP-15C for maximum possible precision. By short I mean no more than 30 or 40 RPN steps.

Late to the party and a bit surprised nobody came up with the following solution:

HP-15C

Code:

001 - 42,21,15  LBL E      
002 -       15  1/x        
003 -    43 32  RTN

1
ENTER
2
∫ E

With FIX 9 it takes about a minute to get:

0.693147181

While CLEAR PREFIX will show:

6931471805

We can compare it to the correct value:

2
LN

0.693147181

But CLEAR PREFIX shows:

6931471806


We can repeat the calculation with FIX 4 and get after about 15 seconds:

0.6931

Here CLEAR PREFIX shows:

6931456394

Interestingly the result is the same with FIX 3 and FIX 5.


HP-42S

Code:

00 { 14-Byte Prgm }
01▸LBL "LN2"
02 MVAR "x"
03 RCL "x"
04 1/X
05 END

We can use ∫ f(x) with the following parameters:

LLIM = 1
ULIM = 2
ACC = 1E-34

The result on Free42 after about 1 second is:

6.931471805599453094172321214581642e-1

Compare it to the correct value:

6.931471805599453094172321214581766e-1

---

While these solutions may not be considered fast enough they are probably the shortest.

(10-12-2019 12:32 AM)Gerson W. Barbosa Wrote:  without using any built-in function (LN, LOG, ATANH...)

I'm aware that calculating an integral is a built-in function, but then calculating the square root is as well.
Thus I assumed that this rule should prevent a trivial way to do it.
But then you may consider this a trivial solution as well.

[Gerson W. Barbosa]

(06-12-2022 12:17 PM)Thomas Klemm Wrote:  

(10-12-2019 12:32 AM)Gerson W. Barbosa Wrote:  By fast I mean under 10 seconds on the HP-15C for maximum possible precision. By short I mean no more than 30 or 40 RPN steps.

Late to the party and a bit surprised nobody came up with the following solution:

HP-15C

Code:

001 - 42,21,15  LBL E      
002 -       15  1/x        
003 -    43 32  RTN

1
ENTER
2
∫ E

With FIX 9 it takes about a minute to get:

0.693147181

While CLEAR PREFIX will show:

6931471805

We can compare it to the correct value:

2
LN

0.693147181

One minute is way longer than 10 seconds, but not too long. Also, on my HP-15C LE it would take no more than half a second.

(06-12-2022 12:17 PM)Thomas Klemm Wrote:  

(10-12-2019 12:32 AM)Gerson W. Barbosa Wrote:  without using any built-in function (LN, LOG, ATANH...)

I'm aware that calculating an integral is a built-in function, but then calculating the square root is as well.
Thus I assumed that this rule should prevent a trivial way to do it.
But then you may consider this a trivial solution as well.

I am not good at setting up rules. Loopholes are always likely to occur. Too good I left Law School after just one week :-)

Now that we are back at this topic, early in the beginning of the pandemic I was asked to find a fast way to compute the binary log of 10 (3.32192809...). I don't know any specific binary log algorithm, but since Albert Chan had already provided a method for instantly computing the natural logarithm of 2, I came up with a solution that yields slightly more than 8 digits per iteration. Probably there are more efficient ways to do it, but I will publish the DECIMAL BASIC code and the first 1000 digits of the constant here just for the record. Too few remarks as I haven't kept any notes.

Code:


OPTION ARITHMETIC DECIMAL_HIGH
LET n = 2^24
LET nd = 1000
PRINT "LB(10) ="
LET t = TIME
LET g = 0.6931471805599453                            !'ln(2) 16 digits
LET h = g/n
LET a = h*h
LET x = 1
LET s = 1 - h                                         !'s = exp(-h)
FOR k = 3 TO 111 STEP 2                               !'loops = 55
   LET x = x * a/(k*k-k)
   LET s = s + x*(k-h)                                !'2 terms at a time
NEXT k
LET s = s^n - 0.5                                     !'s = exp(-g) - 0.5
LET x = s/(s+1)                                       !'ln(1+2s) = 2*atanh(x)
LET a = x*x
LET s = x
FOR k = 3 TO 57 STEP 2                                !'loops = 28
   LET x = x*a
   LET s = s + x/k
NEXT k
LET ln = g + 2*s                                      !ln(2) = g + 2*atanh(x)
LET p1 = 1E643
LET p2 = 2^2136
LET xx = (p2 - p1)/(p2 + p1)
LET x2 = xx*xx
LET sm = 0
LET ni = nd/(LOG(1/xx)/LOG(10)) - 1/2 
FOR nn = 1 TO ni STEP 2
   LET sm = sm + xx/nn
   LET xx = xx*x2
NEXT nn
LET lb = (2136 - (sm + sm)/ln)/643                    !LB(10) = (2136 - 2*sm/ln(2))/643
LET r = TIME - t
LET r$ = STR$(lb - INT(lb))
PRINT
PRINT STR$(INT(lb));
PRINT r$(0:1);
FOR i = 2 TO nd + 1
   PRINT r$(i:i);
   IF MOD((i - 1),10) = 0 THEN PRINT " ";
   IF MOD((i - 1),50) = 0 THEN 
      PRINT
      PRINT REPEAT$(" ",1 + LEN(STR$(INT(lb))));
   END IF
NEXT i
IF MOD (i - 2,50) <> 0  OR nd = 0 THEN PRINT
PRINT 
PRINT "Runtime: ";
PRINT  USING "0.##": r;
PRINT " seconds"
END

--------------------------------------------------------

LB(10) =

3.3219280948 8736234787 0319429489 3901758648 3139302458 
  0612054756 3958159347 7660862521 5850139743 3593701550 
  9965737171 0250251826 8240969842 6352688827 5302772998 
  6553938519 5135265750 5568643017 6091900248 9166694143 
  3374011903 1241873751 0971586646 7540179189 6558067358 
  3077968843 2725883274 9925224489 0238355997 6417394137 
  9280097727 5668635547 7901486745 0578458847 8027104225 
  4560972234 6579569554 1537019157 6411717792 4716513500 
  2392112714 7339361440 7233972115 7485100709 4987891658 
  8808313221 9480679329 8232325931 1950671399 5078370033 
  6734248070 6635275008 4069176263 8625354688 0153686216 
  1841886085 8994835381 3214998930 2704417920 7865922601 
  8229653715 7536723966 0695116486 8368466238 5850848606 
  2990542699 4692791162 7320613400 6446704847 6340704373 
  5233674221 2830896703 6457909216 7721909021 4219621424 
  5744465852 4535948448 8154834592 5142954093 7353906549 
  4486327792 9842429159 1181131163 2981257694 5019815750 
  3792185538 4878203551 6019737827 7288881759 8743328660 
  7271239382 5202213332 8052551248 8274344488 4245316546 
  5061241489 1822867932 5266429281 1659922851 627345082 
  
Runtime: 0.03 seconds

--------------------------------------------------------

Best regards,

Gerson.

[Thomas Klemm]

This program for the HP-42S uses formula (19):

\(
\begin{align}
\ln 2= \sum_{k=1}^{\infty} \frac{1}{k \cdot 2^k}
\end{align}
\)

Code:

00 { 14-Byte Prgm }
01 0
02▸LBL 00
03 RCL ST Y
04 1/X
05 +
06 2
07 ÷
08 DSE ST Y
09 GTO 00
10 END

37
R/S

0.693147180559758517350355291780403
0.6931471805599453094172321214581766


105
R/S

0.6931471805599453094172321214581765
0.6931471805599453094172321214581766


I have no idea how long it takes on a real HP-42S.

[C.Ret]

(06-19-2022 02:10 PM)Thomas Klemm Wrote:  This program for the HP-42S uses formula (19):

\(
\begin{align}
\ln 2= \sum_{k=1}^{\infty} \frac{1}{k \cdot 2^k}
\end{align}
\)

Code:

.. HP42S code ...

I have no idea how long it takes on a real HP-42S.

I like the way Thomas get rid of the \( 2^k \) factor by using a factorized division by 2 in his code. This spare several program steps and globally speed up all the process.

So, I immediately convert his program for my HP-15C:

Code:

001 - {    44,25 }    STO I
002 - {        0 }     0
003 - { 42,21, 0 }  f LBL 0
002 - {    45 25 }    RCL I
005 - {       33 }    1/x
006 - {       40 }     +
007 - {        2 }     2
008 - {       10 }     ÷
009 - { 42, 5,25 }  f DSE I
010 - {     22 0 }    GTO 0

The fastest and accurate result obtained with 27 as argument :
f FIX 9 f PRGM 27 R/S display 0.693147181 in least than 23"97.
The f PREFIX reveals digits 6931471805

[Thomas Klemm]

We can also use Paul's idea and use the Taylor series of:

\(
\begin{align}
\log\left(\frac{1+x}{1-x}\right)
&= 2 x + 2 x^3 \frac{1}{3} + 2 x^5 \frac{1}{5} + 2 x^7 \frac{1}{7} + 2 x^9 \frac{1}{9} + 2 x^{11} \frac{1}{11} + \cdots \\
&= 2 x \left[1 + x^2 \frac{1}{3} + x^4 \frac{1}{5} + x^6 \frac{1}{7} + x^8 \frac{1}{9} + x^{10} \frac{1}{11} + \cdots\right] \\
\end{align}
\)

For \(x\) we solve the equation:

\(
\begin{align}
\frac{1+x}{1-x} = \sqrt[16]{2}
\end{align}
\)

This leads to:

\(
\begin{align}
x = \frac{\sqrt[16]{2} - 1}{\sqrt[16]{2} + 1}
\end{align}
\)

We could maybe increase \(16\) to further reduce \(x\) and thus speed up the summation.
However we experience cancelation in the numerator.
I haven't checked the influence yet.
Thus for now \(16\) is just a wild guess.


Initialization

R00: \(x^2\)
R01: \(32x\)

1
ENTER
2
SQRT
SQRT
SQRT
SQRT
+
LASTX
1
-
X<>Y
/

X↑2
STO 00

LASTX
32
*
STO 01


Programs

HP-42S:

Code:

00 { 18-Byte Prgm }
01 0
02▸LBL 00
03 RCL× 00
04 RCL ST Y
05 IP
06 1/X
07 +
08 DSE ST Y
09 GTO 00
10 RCL× 01
11 END

HP-15C:

Code:

   000 {             } 
   001 {       44 25 } STO I
   002 {       43 35 } g CLx
   003 {    42 21  0 } f LBL 0
   004 {    45 20  0 } RCL × 0
   005 {       45 25 } RCL I
   006 {       43 44 } g INT
   007 {          15 } 1/x
   008 {          40 } +
   009 {    42  5 25 } f DSE I
   010 {       22  0 } GTO 0
   011 {    45 20  1 } RCL × 1

Examples

3.00002
R/S

0.6931471500554851572008806023721818


5.00002
R/S

0.693147180549725015633085213078373


7.00002
R/S

0.6931471805599415808474258676300801


9.00002
R/S

0.6931471805599453079863200827178153


11.00002
R/S

0.6931471805599453094166642087519808


13.00002
R/S

0.6931471805599453094172318905964119


15.00002
R/S

0.6931471805599453094172321213626303


17.00002
R/S

0.6931471805599453094172321214581366


19.00002
R/S

0.6931471805599453094172321214581768


2
LN

0.6931471805599453094172321214581766

[C.Ret]

(10-12-2019 01:30 PM)ttw Wrote:  Perhaps the methods used here may be of interest.

http://elib.mi.sanu.ac.rs/files/journals...tm1212.pdf

In this publication, the five point Gauss-Legendre Log (5P-GLLOG) lead to this approximation :

\( ln(1+x)\approx \frac{7560x+15120x^2+9870x^3+2310x^4+137x^5}{7560+18900x+16800x^2+6300x^3+900x^4​+30x^5} \)

One may deduce that
\( ln(1+1)\approx \frac{7560+15120+9870+2310+137}{7560+18900+16800+6300+900+30} \)
\( ln(2)\approx \frac{34977}{50490}=0.\underline{6931471}57853 \)

A better fraction may be proposed only base of the numeric value without any further mathematical or fundamental consideration:
\( ln(2)\approx \frac{261740}{377611}= 0.\underline{6931471805}641255154113624867919631578529227167640773176628858799 \cdots \)
Best numeric approximation I get on my HP-71B using FRAC$(LOG(2))
which coincidentally correspond to
\(ln(2)\approx[0;1,2,3,1,6,3,1,1,2,1,1,1,1,3,10]\)

[Albert Chan]

Based from [VA]SRC012C

(12-01-2022 08:27 AM)Werner Wrote:  Using the original asymptotic series H(n) = ln(n) + gamma + 1/(2*n) - 1/(12*n^2) + 1/(120*n^4) - ...
results in: H(n-1)-H(n/2-1) - ln(2) = 1/(2n) + 1/(4n^2) - 1/(8n^4) + 1/(4n^6) - 17/(16n^8) + ..

Sum a shell of b-bits integer reciprocals, subtract RHS, we have estimate for ln(2)
Example: for 2-bits integer, 2 and 3:

1/2 + 1/3 - 1/8 - 1/64 + 1/2048 - 1/16384 + 17/1048576 = 2180467/3145728 ≈ 0.693151791890462

Add 1 more asymptote term for good measure, we have:

CAS> LN2(b) := sum(1/k,k,2^(b-1),2^b-1) - horner([7936,-272,16,-2,1,1/2^(b+1)],1/2^(2b+2))

CAS> LN2(2)

8721775/12582912 ≈ 0.69314440091451

CAS> LN2(3)

312589669081/450971566080 ≈ 0.693147179539803

CAS> LN2(4)

27463878294130001/39622001018535936 ≈ 0.693147180559656

Double the bits, we doubled number of correct digits.
OTTH, double the bits required squaring of number of terms, so this is not efficient.

Update: above off by factor of 2, but conclusion still hold.

2^(2b-1) = (2^(b-1))^2 * 2