Threaded Mode | Linear Mode

Namir · 03-10-2014, 05:48 PM

finding the minimmum/maximum of a function using Newton's method.

Memory Map

R0 = tolerance
R1 = x
R2 = h
R3 = f(x), f''(x)
R4 = f(x+h)
R5 = f(x-h)

HP-67 Implementation

Code:

1 LBL A

2 STO 1     # store guess for optimum x  

3 X<>Y

4 STO 0     # store tolerance

5 LBL 0

6 RCL 1     # display current guess for optimum x

7 PAUSE

8 ABS

9 1

10 +

11 EEX

12 CHS

13 3

14 *

15 STO 2     # calculate and store h

16 RCL 1

17 GSB E

18 STO 3     # Calculate and store f(x)

19 RCL 1

20 RCL 2

21 +

22 GSB E

23 STO 4     # Calculate and store f(x+h)

24 RCL 1

25 RCL 2

26 -

27 GSB E

28 STO 5     # Calculate and store f(x-h)

29 RCL 4

30 +

31 RCL 3

32 2

33 *

34 -

35 RCL 2

36 X^2

37 /

38 STO 3     # calculate and store f''(x)

39 RCL 4

40 RCL 5

41 -

42 RCL 2

43 /

44 2         

45 /         # calculate f'(x)

46 RCL 3

47 /         # calculate diff

48 STO- 1    # x = x - diff

49 ABS

50 RCL 0

51 X<=Y?

52 GTO 0

53 RCL 1

54 RTN

55 LBL E

56 EXP

57 LASTX

58 X^2

59 3

60 *

61 -

62 RTN

Usage

LBL E has the code for the function whose minimum or maximum you seek.

1) Enter the tolerance and press [ENTER].
2) Enter a guess for the minimum/maximum and pres [A].
3) The program pauses showing intermediate values for the refined min/max values. When the program has reached the tolerance level it displays the value of the X for the min/max value.

Matt Agajanian · 04-26-2014, 12:30 AM

Very nice! Thanks.

Willy R. Kunz · 05-26-2014, 12:14 PM

Thank you for sharing.

A pity the listings don't come with key codes. This would allow me to copy/paste the programs right into my HP-67/97 emulators for the iPad...

Thomas Klemm · (This post was last modified: 05-26-2014 10:31 PM by Thomas Klemm.)

(05-26-2014 12:14 PM)Willy R. Kunz Wrote: A pity the listings don't come with key codes.

Code:

001  31 25 11    LBL A

002     33 01    STO 1       # store guess for optimum x

003     35 52    X<>Y

004     33 00    STO 0       # store tolerance

005  31 25 00    LBL 0

006     34 01    RCL 1       # display current guess for optimum x

007     35 72    PAUSE

008     35 64    ABS

009        01    1

010        61    +

011        43    EEX

012        42    CHS

013        03    3

014        71    *

015     33 02    STO 2       # calculate and store h

016     34 01    RCL 1

017  31 22 15    GSB E

018     33 03    STO 3       # Calculate and store f(x)

019     34 01    RCL 1

020     34 02    RCL 2

021        61    +

022  31 22 15    GSB E

023     33 04    STO 4       # Calculate and store f(x+h)

024     34 01    RCL 1

025     34 02    RCL 2

026        51    -

027  31 22 15    GSB E

028     33 05    STO 5       # Calculate and store f(x-h)

029     34 04    RCL 4

030        61    +

031     34 03    RCL 3

032        02    2

033        71    *

034        51    -

035     34 02    RCL 2

036     32 54    X^2

037        81    /

038     33 03    STO 3       # calculate and store f''(x)

039     34 04    RCL 4

040     34 05    RCL 5

041        51    -

042     34 02    RCL 2

043        81    /

044        02    2

045        81    /           # calculate f'(x)

046     34 03    RCL 3

047        81    /           # calculate diff

048  33 51 01    STO- 1      # x = x - diff

049     35 64    ABS

050     34 00    RCL 0

051     32 71    X<=Y?

052     22 00    GTO 0

053     34 01    RCL 1

054     35 22    RTN

055  31 25 15    LBL E

056     32 52    EXP

057     35 82    LASTX

058     32 54    X^2

059        03    3

060        71    *

061        51    -

062     35 22    RTN

Cheers
Thomas

**rprosperi** · 05-26-2014, 10:13 PM

Namir - Nice program

Any reason to take out and use my HP-67 is a good one.

Thanks for that!

Willy R. Kunz · (This post was last modified: 07-11-2014 01:24 PM by Willy R. Kunz.)

(05-26-2014 10:13 PM)Thomas Klemm Wrote:

Code:

001 31 25 11 LBL A 002 33 01 STO 1 # store guess for optimum x 003 35 52 X<>Y 004 33 00 STO 0 # store tolerance ...

Cheers
Thomas

Thanks, Thomas. Very kind of you.

This also prompted me to have a hard look at RPN-67/97's import capabilities. The result being that in the next update the emulators will accept listings without key codes. In fact, you can now enter and edit programs by simply typing calculator function names in the listing area.

All of Namir's programs shown here run without any editing after pasting.

UPDATE:
RPN-67 Pro v2.5 and RPN-97 Pro v1.5 are now available, featuring the new import capabilities.
Also in these updates: for the first time ever, real breakpoints to help you debug your program.

PedroLeiva · 05-21-2016, 10:50 PM

(05-26-2014 10:13 PM)Thomas Klemm Wrote: [quote='Willy R. Kunz' pid='12139' dateline='1401106474']
A pity the listings don't come with key codes.

Dear Thomas:
I am using HP 67 SD, I got for Tolerance=5x10-4 and a guessing initial number of 9, a final figure of 2,833144107. Is this correct ?
Now some basic questions please:
1- What equation represents steps 55-62? How I entry an equation to be evaluated?
2- From my example 2,833 is a Max. or a min.? How to obtain both if corresponds?

Thank you in advance, Pedro

Dieter · 05-22-2016, 01:11 AM

(05-21-2016 10:50 PM)PedroLeiva Wrote: 1- What equation represents steps 55-62? How I entry an equation to be evaluated?

As usual, you enter f(x) at LBL E. The argument x is expected in the X-register. Here f(x) obviously is e^x–3x².

(05-21-2016 10:50 PM)PedroLeiva Wrote: 2- From my example 2,833 is a Max. or a min.? How to obtain both if corresponds?

Since f"(x) is stored in R3, you could simply inpect this value (RCL 3) an see whether it's positive (minimum) or negative (maximum).

Here the function has a minimum at x=2,833147892... (solution of f'(x) = e^x–6x = 0).
At ths point f"(x) = e^x–6 is 10,998887... > 0, so this is a minimum. The content of R3 should be close to the latter value.

Dieter

PedroLeiva · 05-22-2016, 01:26 AM

(05-22-2016 01:11 AM)Dieter Wrote:
(05-21-2016 10:50 PM)PedroLeiva Wrote: 1- What equation represents steps 55-62? How I entry an equation to be evaluated?

As usual, you enter f(x) at LBL E. The argument x is expected in the X-register. Here f(x) obviously is e^x–3x².

(05-21-2016 10:50 PM)PedroLeiva Wrote: 2- From my example 2,833 is a Max. or a min.? How to obtain both if corresponds?

Since f"(x) is stored in R3, you could simply inpect this value (RCL 3) an see whether it's positive (minimum) or negative (maximum).

Here the function has a minimum at x=2,833147892... (solution of f'(x) = e^x–6x = 0).
At ths point f"(x) = e^x–6 is 10,998887... > 0, so this is a minimum. The content of R3 should be close to the latter value.

Dieter

Dieter
It's more complicated than I imagined. Luckily I made the query. Thank you very much, Peter

Dieter · (This post was last modified: 05-22-2016 07:36 PM by Dieter.)

(05-22-2016 01:26 AM)PedroLeiva Wrote: It's more complicated than I imagined. Luckily I made the query.

Hmmm... honestly, I do not quite understand the problem.

You simply enter f(x) as a short routine just as you do with any other "classic" HP, be it the 65, the 41, the 67/97 or most others. The argument x can be expected in the X-register, so you simply type what you'd also do in a manual calculation. In this example f(x) is e^x–3x², so it's [e^x] [LstX] [x²] [3] [x] [–]. This is what you enter at LBL E so that f(x) can be calculated by a simple [E] or GSB E.

Checking for a minimum or maximum is easily done by RCL 3 after the program has finished. If the result is positive it's a minimum, and if it's negative there's a maximum. A value of (or very close to) zero would indicate a possible inflection point, but this should not happen as the program divides by R3 so that it would have stopped with an error before.

You could add these lines:

Code:

...  ...

053  RCL 3

054  ENTER

055  ABS

056  /

057  CHS

058  RCL 1

059  RTN

This calculates –sign(f"(x)). So after the program has finished, a simple X<>Y will show –1 (minimum) or +1 (maximum).

Dieter

PedroLeiva · 05-22-2016, 10:35 PM

This calculates –sign(f"(x)). So after the program has finished, a simple X<>Y will show –1 (minimum) or +1 (maximum).
Dieter
[/quote]
The point is that in the original post of this program was no explanation about how to determine Max. or min. (or I did not realize).

From the two options for cheking, I prefer: x<>y (+1, -1), for Max. & min., respectively. This one is more intuitive (+ is related to Max.), the other procedure is inverse.
Thank you for clarifying this math issue, and modifying the program

Pedro

Dieter · 05-23-2016, 01:08 PM

(05-22-2016 10:35 PM)PedroLeiva Wrote: The point is that in the original post of this program was no explanation about how to determine Max. or min. (or I did not realize).

The program implements a simple Newton iteration to find the zeroes of the derivative f'(x). In a minimum or maximum f'(x) must be zero, and this is what the program finds. It evaluates x_new = x – f'(x) / f"(x).

(05-22-2016 10:35 PM)PedroLeiva Wrote: From the two options for cheking, I prefer: x<>y (+1, -1), for Max. & min., respectively.

This is what the proposed code does. It returns +1 for a maximum and –1 for a minimum.

Dieter

PedroLeiva · 05-23-2016, 07:37 PM

Examples:
Equation (began step 061) ex–3x²
Tolerance= 1x10-5
Initial guess: 9……………X= 2,83314411………..min…………..Y= -7,08129358
Initial guess: 0,1…………X= 0,20448151……….Max…………..Y= 1,10145071
For X=0……………………Y= 1
For Y=0……………………X= 0,91

For the original equation, input TOLERANCE [ENTER] 9 [A] you get 2,8331, [x<>y] -1 so you know that is a minimun. Later, [x<>y] [E] you get Y=-7.081 wich is the Y value for the min.
When X=0 [E] you get Y=1, that is to say intercept
The value of X when Y=0 was solved by iteraction

Follow the same procedure for the Maximun
Pedro

bshoring · 05-23-2016, 08:00 PM

(05-26-2014 10:13 PM)Thomas Klemm Wrote:
(05-26-2014 12:14 PM)Willy R. Kunz Wrote: A pity the listings don't come with key codes.

Code:

001 31 25 11 LBL A 002 33 01 STO 1 # store guess for optimum x 003 35 52 X<>Y 004 33 00 STO 0 # store tolerance 005 31 25 00 LBL 0 006 34 01 RCL 1 # display current guess for optimum x 007 35 72 PAUSE 008 35 64 ABS 009 01 1 010 61 + 011 43 EEX 012 42 CHS 013 03 3 014 71 * 015 33 02 STO 2 # calculate and store h 016 34 01 RCL 1 017 31 22 15 GSB E 018 33 03 STO 3 # Calculate and store f(x) 019 34 01 RCL 1 020 34 02 RCL 2 021 61 + 022 31 22 15 GSB E 023 33 04 STO 4 # Calculate and store f(x+h) 024 34 01 RCL 1 025 34 02 RCL 2 026 51 - 027 31 22 15 GSB E 028 33 05 STO 5 # Calculate and store f(x-h) 029 34 04 RCL 4 030 61 + 031 34 03 RCL 3 032 02 2 033 71 * 034 51 - 035 34 02 RCL 2 036 32 54 X^2 037 81 / 038 33 03 STO 3 # calculate and store f''(x) 039 34 04 RCL 4 040 34 05 RCL 5 041 51 - 042 34 02 RCL 2 043 81 / 044 02 2 045 81 / # calculate f'(x) 046 34 03 RCL 3 047 81 / # calculate diff 048 33 51 01 STO- 1 # x = x - diff 049 35 64 ABS 050 34 00 RCL 0 051 32 71 X<=Y? 052 22 00 GTO 0 053 34 01 RCL 1 054 35 22 RTN 055 31 25 15 LBL E 056 32 52 EXP 057 35 82 LASTX 058 32 54 X^2 059 03 3 060 71 * 061 51 - 062 35 22 RTN

Cheers
Thomas

Thanks Thomas! I was able to copy and paste to my iPhone's RPN-67 SD and it works like a charm!

bshoring · 05-24-2016, 06:47 AM

I found a pretty good article about the concept at
http://www.themathpage.com/acalc/max.htm

It helped me to understand what's happening with this program.

bshoring · 05-28-2016, 08:09 PM

I found a pretty good article about the concept at
http://www.themathpage.com/acalc/max.htm

It helped me to understand what's happening with this program.

PedroLeiva · 05-28-2016, 09:43 PM

(05-28-2016 08:09 PM)bshoring Wrote: I found a pretty good article about the concept at
http://www.themathpage.com/acalc/max.htm

It helped me to understand what's happening with this program.

Excellent webpage, also for the others 14th. lessons. Thank´s for sharing
Pedro