# HP Forums

Full Version: Small RPN exponential routine
You're currently viewing a stripped down version of our content. View the full version with proper formatting.
Pages: 1 2
This would make more sense on the HP-16C, which lacks exponentiation, but the equivalent program is 23 steps long.
I wrote a 42S version as an optimization exercise and managed to save five steps and two registers (but only two bytes), when compared to the first version. Perhaps more optimization is still possible, but I'll stop here.
I'll publish the HP-16C version, in case someone is interested (on the 16C this should be a first step towards a usable exponential function implementation).

Code:
 00 { 29-Byte Prgm } 01>LBL "EX" 02 RCL ST Y 03 1 04>LBL 00 05 RCL× ST Z 06 STO+ ST Y 07 X<>Y 08 RCL× ST T 09 X<>Y 10 DSE ST Z 11 GTO 00 12 ÷ 13 RCL÷ ST Z 14 .END.

Examples:

1 ENTER 14 XEQ EX --> 2.71828182846

3.14159265359 ENTER 22 --> 23.1406926328

1000 ENTER 1400 --> 1.97007111402e434 ( Free 42 only! )
(01-02-2015 10:25 PM)Gerson W. Barbosa Wrote: [ -> ]Perhaps more optimization is still possible, but I'll stop here.

Is there a reason for not using the obvious?
Code:
00 { 22-Byte Prgm } 01>LBL "EX" 02 1 03>LBL 00 04 RCL× ST Z 05 RCL÷ ST Y 06 1 07 + 08 DSE ST Y 09 GTO 00 10 END

Cheers
Thomas
(01-03-2015 08:42 AM)Thomas Klemm Wrote: [ -> ]Is there a reason for not using the obvious?

Yes, dumbness...

D'oh!

Cheers,

Gerson.
This is my first attempt at natural logarithm using

ln(1 + x) = x - x^2/2 + x^3/3 - x^4/4 + x^5/5 -+ ..., x < 1

There surely is a more obvious method which I cannot see, so I won't bother optimizing it. I am primarily interested in implementing these functions on TurboBCD, but rather than trying to reinvent the wheel I'll look for a ready-made solution.

Code:
 00 { 43-Byte Prgm } 01>LBL "LX" 02 NOT 03 +/- 04 LASTX 05 1 06 RCL- ST T 07 RCL÷ ST Z 08 LASTX 09 1 10>LBL 00 11 RCL÷ ST T 12 RCL+ ST Z 13 RCL× ST Y 14 STO ST Z 15 Rv 16 1 17 DSE ST T 18 GTO 00 19 RCL ST Z 20 +/- 21 .END.

0.52 ENTER 35 --> -6.53926467407e-1

There are two problems, though:

1) The series does not converge for arguments greater than 1;

2) More and more terms are required as arguments approach to zero. For instance, ln(0.052) will require 525 terms for 12-digit accuracy.

In these cases, conversions to and from decimal logarithms will help:

ln(0.52) = -6.53926467407e-1

2.30258509299 / --> -2.83996656366e-1

1 - 2.30258509299 * --> -2.95651156041 = ln(0.052)

ln(0.52) = -6.53926467407e-1

2.30258509299 / --> -2.83996656366e-1

1 + 2.30258509299 * --> 1.64865862558 = ln(5.2)

ln(0.52) = -6.53926467407e-1

2.30258509299 / --> -2.83996656366e-1

2 + 2.30258509299 * --> 3.95124371856 = ln(52)

ln(0.52) = -6.53926467407e-1

2.30258509299 / --> -2.83996656366e-1

3 + 2.30258509299 * --> 6.25382881155 = ln(520)
(01-03-2015 05:31 PM)Gerson W. Barbosa Wrote: [ -> ]There surely is a more obvious method which I cannot see, so I won't bother optimizing it.

Using the same method but a little shorter:
Code:
00 { 24-Byte Prgm } 01>LBL "LX" 02 1 03 STO- ST Z 04 CLX 05>LBL 00 06 RCL Y 07 1/X 08 X<>Y 09 - 10 RCL× ST Z 11 DSE ST Y 12 GTO 00 13 END

Cheers
Thomas
(01-03-2015 08:33 PM)Thomas Klemm Wrote: [ -> ]
(01-03-2015 05:31 PM)Gerson W. Barbosa Wrote: [ -> ]There surely is a more obvious method which I cannot see, so I won't bother optimizing it.

Using the same method but a little shorter:
Code:
00 { 24-Byte Prgm } 01>LBL "LX" 02 1 03 STO- ST Z 04 CLX 05>LBL 00 06 RCL Y 07 1/X 08 X<>Y 09 - 10 RCL× ST Z 11 DSE ST Y 12 GTO 00 13 END

Cheers
Thomas

A "little shorter"? I appreciate your sense of humor :-)

Definitely I should have made a second attempt in the first place. I've rewritten it looking at the expression with five terms of my example and obtained the same number of steps, but four bytes longer. Despite having one extra step inside the loop, yours appears to be a little faster (only one measurement before I got "Memory Clear" due to flat batteries).

Cheers,

Gerson.

ln(0.52) = -(0.48*(1 + 0.48*(1/2 + 0.48*(1/3 + 0.48*(1/4 + 0.48/5)))))

Code:
 00 { 28-Byte Prgm } 01>LBL "LX" 02 1 03 RCL- ST Z 04 0 05>LBL 00 06 1 07 RCL÷ ST T 08 + 09 RCL× ST Y 10 DSE ST Z 11 GTO 00 12 +/- 13 .END.

P.S.: One step less, but still one byte longer, using one of your ideas:

Code:
 00 { 25-Byte Prgm } 01>LBL "LX" 02 1 03 RCL- ST Z 04 0 05>LBL 00 06 RCL ST Z 07 1/X 08 - 09 RCL× ST Y 10 DSE ST Z 11 GTO 00 12 END
Gerson, Is your exponential routine similar to Toth's? (listed below)

Code:
001 - 43,22, d  LBL D 002 - 44 8      STO 8 003 - 1         1 004 - 44 32     STO I 005 - 44 9      STO 9 006 - 20        × 007 - 43,22, 4  LBL 4 008 - 36        ENTER 009 - 36        ENTER 010 - 45 9      RCL 9 011 - 40        + 012 - 45 9      RCL 9 013 - 43 49     x=y 014 - 43 21     RTN 015 - 33        Rv 016 - 44 9      STO 9 017 - 43 24     ISZ 018 - 33        Rv 019 - 45 8      RCL 8 020 - 20        × 021 - 45 32     RCL I 022 - 10        ÷ 023 - 22 4      GTO 4
It works, but I wonder how. Could someone explain? Thanks.
How about HP-41C versions .... pleeeeeeeease!!

:-)

Namir
(01-20-2015 07:57 PM)James Dunn Wrote: [ -> ]Gerson, Is your exponential routine similar to Toth's? (listed below)

Code:
001 - 43,22, d  LBL D 002 - 44 8      STO 8 003 - 1         1 004 - 44 32     STO I 005 - 44 9      STO 9 006 - 20        × 007 - 43,22, 4  LBL 4 008 - 36        ENTER 009 - 36        ENTER 010 - 45 9      RCL 9 011 - 40        + 012 - 45 9      RCL 9 013 - 43 49     x=y 014 - 43 21     RTN 015 - 33        Rv 016 - 44 9      STO 9 017 - 43 24     ISZ 018 - 33        Rv 019 - 45 8      RCL 8 020 - 20        × 021 - 45 32     RCL I 022 - 10        ÷ 023 - 22 4      GTO 4

My HP-16C version:

Code:
 01 g LBL A 02   STO I 03   1 04   STO 0 05 g CLx 06   x<>y 07   ENTER 08   ENTER 09   ENTER 10 g LBL 0 11   RCL 0 12   RCL I 13   * 14   STO 0 15   +  16   * 17 g DSZ 20   GTO 0 19   RCL 0 20   / 21   x<>y 22   / 23 g RTN

Both have approximately the same size and both use Horner's method. My first HP-16C version of Thomas Klemm's more obvious method is only one step shorter, but I surely haven't converted it properly.

Regards,

Gerson.

P.S.: Contrary to what I stated, Victor Toth's program uses a different method. My bad, sorry!
(01-20-2015 11:48 PM)Gerson W. Barbosa Wrote: [ -> ]My HP-16C version:

Code:
 01 g LBL A 02   STO I 03   1 04   STO 0 05 g CLx 06   x<>y 07   ENTER 08   ENTER 09   ENTER 10 g LBL 0 11   RCL 0 12   RCL I 13   * 14   STO 0 15   +  16   * 17 g DSZ 20   GTO 0 19   RCL 0 20   / 21   x<>y 22   / 23 g RTN

Both have approximately the same size and both use Horner's method. My first HP-16C version of Thomas Klemm's more obvious method is only one step shorter, but I surely haven't converted it properly.

Regards,

Gerson.

I think LBL A needs to be followed by an X<>Y command. I am assuming the value for the exponent is in register Y and the number of iterations is in register X. I tried your code on an HP-41C and it worked after inserting X<>Y following LBL A.

Namir
(01-21-2015 12:34 AM)Namir Wrote: [ -> ]
(01-20-2015 11:48 PM)Gerson W. Barbosa Wrote: [ -> ]My HP-16C version:

Code:
 01 g LBL A 02   STO I 03   1 04   STO 0 05 g CLx 06   x<>y 07   ENTER 08   ENTER 09   ENTER 10 g LBL 0 11   RCL 0 12   RCL I 13   * 14   STO 0 15   +  16   * 17 g DSZ 20   GTO 0 19   RCL 0 20   / 21   x<>y 22   / 23 g RTN

Both have approximately the same size and both use Horner's method. My first HP-16C version of Thomas Klemm's more obvious method is only one step shorter, but I surely haven't converted it properly.

Regards,

Gerson.

I think LBL A needs to be followed by an X<>Y command. I am assuming the value for the exponent is in register Y and the number of iterations is in register X. I tried your code on an HP-41C and it worked after inserting X<>Y following LBL A.

Namir

Namir,

On the 16C, 1 ENTER 12 GSB A returns 2.718281829. No x<>y need. I'll try it on the 41C.

Gerson.

P.S.: Instead of X<>Y after LBL A, try

STO 01 ; assuming that's your index register;
SIGN
STO 00
X<>Y
ENTER
ENTER
ENTER
...

Of course converting Thomas Klemm's program to the HP-41C is better.
(01-20-2015 07:57 PM)James Dunn Wrote: [ -> ]Gerson, Is your exponential routine similar to Toth's?
It works, but I wonder how. Could someone explain? Thanks.

Toth's program calculates the partial sums:

$$1$$
$$1+x$$
$$1+x+\frac{x^2}{2}$$
$$1+x+\frac{x^2}{2}+\frac{x^3}{6}$$
$$1+x+\frac{x^2}{2}+\frac{x^3}{6}+\frac{x^4}{24}$$
(...)

It returns if the value doesn't change anymore. This value is kept in register 9 while the original x is kept in register 8. The terms to be added are calculated recursively from the previous value: $$a_n=\frac{x}{n}a_{n-1}$$ (lines 019-021). The index $$n$$ is kept in register I.

Horner's method let us rewrite the expression:

$$1+x(1+\frac{x}{2}(1+\frac{x}{3}(1+\frac{x}{4}(...)))$$

My program calculates that inside out using classic RPN style. Thus the loop is inverted: the index is going down to 0. From a numerical point of view that's favorable as it reduces cancelation: the values that are added are of a similar size.

Best is to step through the program using special values like 2 that are easy to follow.

HTH
Thomas
(01-20-2015 11:48 PM)Gerson W. Barbosa Wrote: [ -> ]My first HP-16C version of Thomas Klemm's more obvious method is only one step shorter, but I surely haven't converted it properly.

This is how I would write my program for the HP-16C:
Code:
001 - 43,22, A  LBL A 002 - 44 32     STO I 003 - 33        X<>Y 004 - 44 0      STO 0 005 - 1         1 006 - 43,22, 0  LBL 0 007 - 45 0      RCL 0 008 - 20        × 009 - 45 32     RCL I 010 - 10        ÷ 011 - 1         1 012 - 40        + 013 - 43 23     DSZ 014 - 22 0      GTO 0 015 - 43 21     RTN

Cheers
Thomas
(01-21-2015 03:18 AM)Thomas Klemm Wrote: [ -> ]
(01-20-2015 11:48 PM)Gerson W. Barbosa Wrote: [ -> ]My first HP-16C version of Thomas Klemm's more obvious method is only one step shorter, but I surely haven't converted it properly.

This is how I would write my program for the HP-16C:
Code:
001 - 43,22, A  LBL A 002 - 44 32     STO I 003 - 33        X<>Y 004 - 44 0      STO 0 005 - 1         1 006 - 43,22, 0  LBL 0 007 - 45 0      RCL 0 008 - 20        × 009 - 45 32     RCL I 010 - 10        ÷ 011 - 1         1 012 - 40        + 013 - 43 23     DSZ 014 - 22 0      GTO 0 015 - 43 21     RTN

Cheers
Thomas

I'm glad my 22-step version has gone straight to the trash bin :-)

Cheers,

Gerson.
(01-20-2015 09:24 PM)Namir Wrote: [ -> ]How about HP-41C versions .... pleeeeeeeease!!

:-)

Namir

Conversions from Thomas Klemm's HP-42S programs above are straightforward:

Code:
 01 LBL 'LN                  01 LBL 'EX 02 1                        02 1             03 ST- Z                    03 LBL 00 04 CLX                      04 RCL Z 05 LBL 00                   05 * 06 RCL Y                    06 RCL Y 07 1/X                      07 / 08 X<>Y                     08 1 09 -                        09 + 10 RCL Z                    10 DSE Y 11 *                        11 GTO 00 12 DSE Y                    12 END 13 GTO 00 14 END

Gerson.
Thank you for the HP-41C Code. Now for a slightly more difficult request. How about a version fo rthe HP-67?

Namir
(01-22-2015 08:41 AM)Namir Wrote: [ -> ]Thank you for the HP-41C Code. Now for a slightly more difficult request. How about a version fo rthe HP-67?

Namir

No HP-67 here, but I have a printed manual. These should be equivalent to the above HP-41C programs and hopefully work on the HP-67:

Code:
 001 h LBL A                    016 h LBL B 002 h ST I                     017 h ST I 003 h x<>y                     018 h x<>y 004   STO 1                    019   1 005   1                        020   - 006 h LBL 1                    021   STO 1 007   RCL 1                    022   0 008   *                        023 h LBL 2 009 h RC I                     024 h RC I 010   /                        025   1/x 011   1                        026 h x<>y 012   +                        027   - 013 f DSZ                      028   RCL 1 014   GTO 1                    029   * 015 h RTN                      030 f DSZ                                031   GTO 2                                032 h RTN

Gerson.
(01-22-2015 06:59 PM)Gerson W. Barbosa Wrote: [ -> ]
(01-22-2015 08:41 AM)Namir Wrote: [ -> ]Thank you for the HP-41C Code. Now for a slightly more difficult request. How about a version fo rthe HP-67?

Namir

No HP-67 here, but I have a printed manual. These should be equivalent to the above HP-41C programs and hopefully work on the HP-67:

Code:
 001 h LBL A                    016 h LBL B 002 h ST I                     017 h ST I 003 h x<>y                     018 h x<>y 004   STO 1                    019   1 005   1                        020   - 006 h LBL 1                    021   STO 1 007   RCL 1                    022   0 008   *                        023 h LBL 2 009 h RC I                     024 h RC I 010   /                        025   1/x 011   1                        026 h x<>y 012   +                        027   - 013 f DSZ                      028   RCL 1 014   GTO 1                    029   * 015 h RTN                      030 f DSZ                                031   GTO 2                                032 h RTN

Gerson.

Gerson,

Thank you so much for the HP-67 code.

Have you considered attending an HHC conference in the US? It will be a special treat to meet you in person. I am sure you will have a time of your life being among peers who respect you and your contributions here.

Namir
(01-23-2015 05:11 PM)Namir Wrote: [ -> ]Gerson,

Thank you so much for the HP-67 code.

Have you considered attending an HHC conference in the US? It will be a special treat to meet you in person. I am sure you will have a time of your life being among peers who respect you and your contributions here.

Namir

Namir,

Thank you very much for your kind words, despite my ever fading contributions.
Yes, I've already considered attendind one HHC conference. It would be a pleasure meeting you all in person. Perhaps one of these years :-)

Thanks again for the invitation!

Gerson.
Is this the origin for the LPN1 function? (I need to be here more often)
Pages: 1 2
Reference URL's
• HP Forums: https://www.hpmuseum.org/forum/index.php
• :