(35S) Extremum search - Function minimization

[stephane]

We often need to optimize a function. Unfortunately, the equation solver of the HP 35s calculator is not able to find an extremum of a function.

If you can't calculate the derivative of the function, or don't want to calculate it, you can then run the program below, based on the Golden section search algorithm.
The program computes the minimum of a function on a given interval of values [a, b]. The function must be unimodal on this interval, otherwise the algorithm could only find a local minimum, not the global one.

The function to optimize must be stored in the program Y. It must take the x value on the stack (REGX), and it must return the result into the stack too (REGX). The entire stack can be used by the program Y to compute the result of the function.

Example:
y=(x-2)²+3

A program for this function can be:

Code:


Y001 LBL Y
Y002 2
Y003 -
Y004 x²
Y005 3
Y006 +
Y007 RTN

CK=40C4 LN=29

or

Code:


Y001 LBL Y
Y002 STO X
Y003 SQ(X-2)+3
Y004 RTN

CK=2FBC LN=21

The minimization program that follows uses the variables A..F of the calculator, so, don't use them in the equation Y. The other variables are preserved, and are then available.


Usage to minimize a function f(x) stored in Y:

Stack before:
Y: x_min
X: x_max

Run the program: XEQ E (E means Extremum, and it is on the Rdown key with the down arrow showing that we are looking for a minimum)
The program then ask the user for a search tolerance value E:

E?
0.0001

You can change the value if your want, but a lower value will take a longer time of computation.
Press then R/S to begin the minimum search.

Stack after:
Y: x_min'
X: x_max'

x_max' - x_min' is <= given tolerance E.
As x_min' and x_max' are new x bounds for the minimum value of the function, you can run again the program with a smaller value of tolerance E, or you can compute the corresponding value of the function simply keying XEQ Y now.

Code of the program:

Code:


E001 LBL E
E002 SF 10
E003 MINIMIZE Y    CK=4275 LN=10
E004 PSE
E005 CF 10
E006 1E-4
E007 STO E
E008 Rdown
E009 INPUT E
E010 Rdown
E011 STO B
E012 x<>y
E013 STO A
E014 +/-
E015 RCL+ B
E016 0.618
E017 *
E018 +/-
E019 RCL+ B
E020 STO C
E021 XEQ Y001
E022 STO F
E023 RCL B
E024 RCL- A
E025 0.618
E026 *
E027 RCL+ A
E028 STO D
E029 XEQ Y001
E030 RCL F
E031 x<y?
E032 GTO E036
E033 RCL C
E034 RCL B
E035 GTO E038
E036 RCL A
E037 RCL D
E038 REGX
E039 REGZ
E040 -
E041 RCL- E
E042 x>0?
E043 GTO E010
E044 Rdown
E045 RTN

CK=A6AA LN=159

Note: 0.618 comes from the Golden number, and equals 2/(1+sqrt(5))