The Museum of HP Calculators


Complex Functions for the HP-41

This program is Copyright © 2004 by Jean-Marc Baillard and is used here by permission.

This program is supplied without representation or warranty of any kind. Jean-Marc Baillard and The Museum of HP Calculators therefore assume no responsibility and shall have no liability, consequential or otherwise, of any kind arising from the use of this program material or any part thereof.

Overview

-The following programs calculate the sum, difference, product and quotient of 2 complex numbers,
  the square, square root, natural logarithm, natural exponential, inverse of one complex number, trigonometric and hyperbolic functions and their inverses,
  raising a complex number to a real power and raising a complex number to a complex power. ( 23 routines in all )
-Then, 4 root-finding programs are listed.
 

1°) Complex Operations


-The global labels  "Z+"  "Z-"  "Z*"  "Z/"   mean addition, subtraction, multiplication and division.
  "Z^2"  "SQRTZ"  "LNZ"  "E^Z"  "1/Z"   -----  square, square root, logarithm, exponential, reciprocal.
  "SINZ"  "COSZ"  "TANZ"  "ASINZ"  "ACOSZ"  "ATANZ"  correspond to sine, cosine, tangent and their inverses.
  "SHZ"  "CHZ"  "THZ"  "ASHZ"  "ACHZ"  "ATHZ"  correspond to hyperbolic sine, hyperbolic cosine, hyperbolic tangent and their inverses:
-I've used the French notations ( sh z instead of sinh z  ... etc ... ) which have the advantage of saving several bytes and keystrokes.
  "Z^X" and "Z^Z"  stand for raising a complex number to a real power and a complex power respectively.

-These programs use no data register, no synthetic register and no external subroutine.
-However, for instance, "TANZ"  calls  "THZ"  as a subroutine since  tan z = i.tanh (-i.z) and similarly for many other complex functions.
-In these cases, local labels have been inserted after the global labels  ( for example,  XEQ 03  LBL 03 instead of  XEQ "THZ",
  which saves 1 byte and executes faster ) and so on ...
-Furthermore, some series of instructions are executed several times ( like lines 82 to 90 ), also leading to local subroutines.
-Lines 120 & 139 are only useful te restore the DEG mode: they may be deleted if RAD is your favorite angular mode.

-To get the magnitude of a complex number, simply press  R-P.

Data Registers: /
Flags: /
Subroutines:
  /


001  LBL "SINZ"
002  X<>Y
003  CHS
004  XEQ 01
005  CHS
006  X<>Y
007  RTN
008  LBL "TANZ"
009  CHS
010  X<>Y
011  XEQ 03
012  X<>Y
013  CHS
014  RTN
015  LBL "ASINZ"
016  X<>Y
017  CHS
018  XEQ 05
019  CHS
020  X<>Y
021  RTN
022  LBL "ATANZ"
023  CHS
024  X<>Y
025  XEQ 06
026  X<>Y
027  CHS
028  RTN
029  LBL "SHZ"
030  LBL 01
031  XEQ 00
032  XEQ 01
033  GTO 10
034  LBL "COSZ"
035  X<>Y
036  CHS
037  LBL "CHZ"
038  XEQ 00
039  XEQ 02
040  GTO 10
041  LBL "ATHZ"
042  LBL 06
043  RCL X
044  1
045  ST+ Z
046  X<>Y
047  -
048  R^
049  CHS
050  X<>Y
051  XEQ 07
052  XEQ 12
053  LBL 10
054  2
055  ST/ Z
056  /
057  RTN
058  LBL "THZ"
059  LBL 03
060  2
061  ST* Z
062  *
063  XEQ 08
064  RCL X
065  1
066  ST- Z
067  +
068  R^
069  X<>Y
070  LBL "Z/"
071  LBL 07
072  R-P
073  R^
074  R^
075  R-P
076  R^
077  ST- Z
078  X<> T
079  /
080  P-R
081  RTN
082  LBL 00
083  RCL Y
084  RCL Y
085  XEQ 08
086  R^
087  CHS
088  R^
089  CHS
090  GTO 08
091  LBL "ACOSZ"
092  XEQ 04
093  CHS
094  X<>Y
095  GTO 13
096  LBL "ACHZ"
097  LBL 04
098  XEQ 14
099  CHS
100  XEQ 09
101  LBL 13
102  ENTER^
103  SIGN
104  ST* Z
105  *
106  RTN
107  LBL "ASHZ"
108  LBL 05
109  XEQ 14
110  LBL 09
111  ST+ L
112  X<> L
113  XEQ 11
114  XEQ 02
115  LBL "LNZ"
116  LBL 12
117  RAD
118  R-P
119  LN
120  DEG
121  RTN
122  LBL 14
123  RCL Y
124  RCL Y
125  XEQ 05
126  SIGN
127  ST* X
128  RTN
129  LBL "Z^Z"
130  R^
131  R^
132  XEQ 12
133  XEQ 03
134  LBL "E^Z"
135  LBL 08
136  RAD
137  E^X
138  P-R
139  DEG
140  RTN
141  LBL "Z^X"
142  RDN
143  R-P
144  R^
145  ST* Z
146  Y^X
147  P-R
148  RTN
149  LBL "Z-"
150  LBL 01
151  CHS
152  X<>Y
153  CHS
154  X<>Y
155  LBL "Z+"
156  LBL 02
157  ST+ Z
158  X<> T
159  +
160  X<>Y
161  RTN
162  LBL "Z*"
163  LBL 03
164  R-P
165  R^
166  R^
167  R-P
168  X<>Y
169  ST+ T
170  RDN
171  *
172  P-R
173  RTN
174  LBL "1/Z"
175  R-P
176  1/X
177  X<>Y
178  CHS
179  X<>Y
180  P-R
181  RTN
182  LBL "Z^2"
183  LBL 05
184  R-P
185  X^2
186  X<>Y
187  ST+ X
188  X<>Y
189  P-R
190  RTN
191  LBL "SQRTZ"
192  LBL 11
193  R-P
194  SQRT
195  SIGN
196  ST+ X
197  ST/ Y
198  X<> L
199  P-R
200  END

( 420 bytes / SIZE 000 )


1°) First case:  Function of one complex number


      STACK        INPUTS      OUTPUTS
           Y             y             v
           X             x             u

where   z = x+i.y  ;  f(z) = u+i.v

Example:       z = 2+3.i  
     
-Calculate  z2 ;  z1/2  ;  1/z  ;  exp(z)  ;  Ln(z)  ;  sin z  ;  cos z  ;  tan z ;  Asin z  ; Acos z ; Atan z ;  Sinh z ; Cosh z ; Tanh z ; Asinh z ; Acosh z  ; Atanh z

3  ENTER^   2   XEQ "Z^2"        >>>      -5      X<>Y      12        whence    (2+3.i)2   =   -5+12.i
3  ENTER^   2   XEQ "SQRTZ"  >>>   1.6741  X<>Y   0.8960   --------    (2+3.i)1/2 =   1.6741+0.8960.i
3  ENTER^   2   XEQ "1/Z"         >>>   0.1538  X<>Y  -0.2308   --------   1/(2+3.i)  =    0.1538-0.2308.i
3  ENTER^   2   XEQ "E^Z"        >>>  -7.3151  X<>Y   1.0427   --------   exp(2+3.i) =  -7.3151+1.0427.i
3  ENTER^   2   XEQ "LNZ"       >>>   1.2825  X<>Y    0.9828   --------   Ln(2+3.i)  =   1.2825+0.9828.i
3  ENTER^   2   XEQ "SINZ"     >>>    9.1545  X<>Y   -4.1689  --------   Sin(2+3.i)  =   9.1545-4.1689.i
3  ENTER^   2   XEQ "COSZ"    >>>   -4.1896  X<>Y  -9.1092  --------   Cos(2+3.i) =  -4.1896-9.1092.i
3  ENTER^   2   XEQ "TANZ"    >>>   -0.0038  X<>Y   1.0032  --------    Tan(2+3.i) =  -0.0038+1.0032.i
3  ENTER^   2   XEQ "ASINZ"   >>>    0.5707  X<>Y   1.9834  --------   Asin(2+3.i) =   0.5707+1.9834.i
3  ENTER^   2   XEQ "ACOSZ"  >>>   1.0001  X<>Y  -1.9834  --------   Acos(2+3.i) =  1.0001-1.9834.i
3  ENTER^   2   XEQ "ATANZ"  >>>   1.4099  X<>Y   0.2291  --------    Atan(2+3.i) =  1.4099+0.2291.i
3  ENTER^   2   XEQ "SHZ"        >>>  -3.5906  X<>Y  0.5309  --------    Sinh(2+3.i)  =  -3.5906+0.5309.i
3  ENTER^   2   XEQ "CHZ"       >>>   -3.7245  X<>Y  0.5118  --------    Cosh(2+3.i) = -3.7245+0.5118.i
3  ENTER^   2   XEQ "THZ"       >>>    0.9654   X<>Y  -0.0099  -------    Tanh(2+3.i)  =  0.9654-0.0099.i
3  ENTER^   2   XEQ "ASHZ"    >>>    1.9686   X<>Y   0.9647  --------   Asinh(2+3.i) =  1.9686+0.9647.i
3  ENTER^   2   XEQ "ACHZ"   >>>     1.9834   X<>Y  1.0001  --------    Acosh(2+3.i) = 1.9834+1.0001.i
3  ENTER^   2   XEQ "ATHZ"    >>>    0.1469   X<>Y   1.3390  --------    Atanh(2+3.i) =  0.1469+1.3390.i

Note:    "Z^2"  "SQRTZ"  "1/Z"  "E^Z"  and  "LNZ"   save registers Z and T ( this feature is important when they are called as subroutines by the other functions ).


2°) Second case:  raising a complex number to a real power


      STACK        INPUTS      OUTPUTS
           Z             y             /
           Y             x             v
           X             r             u

where   z = x+i.y  ;  zr = u+i.v

Example:   Evaluate  (2+3.i)1/5

3  ENTER^   2  ENTER^   5   1/X   XEQ "Z^X"   >>>  1.2675   X<>Y   0.2524    whence   (2+3.i)1/5 = 1.2675+0.2524.i

Note:   "Z^X"  saves T-register in  registers Z and T.


3°) Third case:  Function of two complex numbers.


      STACK        INPUTS      OUTPUTS
           T
            y
            /
           Z             x             /
           Y             y'             v
           X             x'             u

where   z = x+i.y  ;  z' = x'+i.y' ; f(z;z') = u+i.v

Example:   z = 2+3.i  ;  z' = 4+7.i  Compute  z+z' ; z-z' ; z.z' ; z/z' ; zz'

3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z+"    gives      6       X<>Y       10      whence    z+z' = 6+10.i
3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z-"     ----      -2       X<>Y       -4      --------    z-z' = -2-4.i
3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z*"     ----    -13       X<>Y       26      --------    z.z' = -13+26.i
3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z/"      ----   0.4462  X<>Y  -0.0308   --------    z/z' = 0.4462-0.0308.i
3  ENTER^   2  ENTER^   7  ENTER^   4  XEQ "Z^Z"   ----  0.1638   X<>Y   0.0583    --------    zz' = 0.1638+0.0583.i


2°) Complex Equations:   f(z) = 0


        a) The "Secant" Method

-This program uses an iterative algorithm.
-Starting with 2 distinct initial guesses  z0 = x0+i.y0 ;  z1 = x1+i.y1 ,  it calculates  zk+1 = zk - f(zk).(zk - zk-1) / [ f(zk) - f(zk-1) ]     ( k = 1 ; 2 ; 3 ; ....... )
-The process is repeated until 2 consecutive approximations are equal.

Data Registers:   • R00 = Function name   ( Global label, maximum of 6 characters. This register is to be initialized before executing "ROOTZ" )
                                 R01 = x   R03 = a   R05 = x'
                                 R02 = y   R04 = b   R06 = y'   where  f(z) = f(x+i.y) = a+i.b
Flag: /
Subroutine:  A program which computes  f(z) = a+i.b  assuming  x is in X-register and y is in Y-register upon entry
                       In other words, the stack must be changed from      Y = y           Y = b
                                                                                                     X = x    to    X = a      where z = x+i.y   and   f(z) = f(x+i.y) = a+i.b


01  LBL "ROOTZ"
02  STO 01
03  X<>Y
04  STO 02
05  R^
06  STO 06
07  R^
08  STO 05
09  XEQ IND 00
10  STO 03
11  X<>Y
12  STO 04
13  LBL 01
14  VIEW 01
15  RCL 02
16  RCL 01
17  XEQ IND 00
18  RCL 04
19  RCL 03
20  R^
21  STO 04
22  ST- Z
23  R^
24  STO 03
25  ST- Z
26  R-P
27  R^
28  R^
29  R-P
30  X<>Y
31  ST- T
32  RDN
33  X#0?
34  /
35  RCL 02
36  ENTER^
37  X<> 06
38  -
39  RCL 01
40  STO L
41  X<> 05
42  ST- L
43  X<> L
44  R-P
45  X<>Y
46  ST+ T
47  RDN
48  *
49  P-R
50  ST+ 01
51  X<>Y
52  ST+ 02
53   E-8
54  LASTX
55  X>Y?
56  GTO 01
57  RCL 02
58  RCL 01
59  END

( 86 bytes / SIZE 007 )


      STACK        INPUTS      OUTPUTS
           T             y0             /
           Z             x0             /
           Y             y1             y
           X             x1             x

Example:   Find a root of      sinh z + z2 + Pi = 0.      First, we load a routine to compute f(z) , for instance:

  LBL "T"
  STO 07
  X<>Y
  STO 08
  X<>Y
  XEQ "SHZ"
  RCL 08
  RCL 07
  XEQ "Z^2"
  XEQ "Z+"
  PI
  +
  RTN

And if we choose   z0 = 0  ;  z1 = 1+i

   T  ASTO 00
   0  ENTER^   ENTER^
   1  ENTER^   XEQ "ROOTZ"  the succesive x-values are displayed and we get   -0.278189856                  ( execution time = 77 seconds )
                                                                                                               X<>Y      1.812880366

Whence     -0.278189856 + 1.812880366.i    is a solution.   ( There are many other solutions like  4.419361546 + 4.697881722.i  ... etc ... )

Note:  -The iterations stop when  zk+1 =  zk  or when   f(zk+1) =  f(zk)  Therefore, check f(z) in registers R03 & R04
 

        b) Quadratic Interpolation

-Starting with 3 initial guesses:    z0 = x0+i.y0 ;  z1 = x1+i.y1 ;  z2 = x2+i.y2 , this program calculates

   A = ( z1-z0 ).f(z2) + ( z0-z2 ).f(z1) + ( z2-z1 ).f(z0)
   B = ( z1-z0 ).( 2.z2-z1-z0 ).f(z2) - ( z0-z2 )2.f(z1) + ( z2-z1 )2.f(z0)
   C = ( z2-z1 ).( z1-z0 ).( z2-z0 ).f(z2)

and   z3 = z2 - 1 / [ B/(2.C) + E.( B2/(4C2) - A/C )1/2 ]      where   E = +1 or -1   ( the sign is choosen to produce the larger denominator )

-The process is repeated until  | zk+1-zk |  is smaller than  10-8  or  C = 0   ( check  f(z) in  R03 & R04 )

Data Registers:   • R00 = Function name   ( Global label, maximum of 6 characters. This register is to be initialized before executing "RTZ2" )
                                 R01 = x   R03 = a   R05 thru R20:  temp
                                 R02 = y   R04 = b   where  f(z) = f(x+i.y) = a+i.b
Flag: /
Subroutines:  "Z+"  "Z*"  "Z^2"  "SQRTZ"
                       and a program which computes  f(z) = a+i.b  assuming  x is in X-register and y is in Y-register upon entry
                       The stack must be changed from      Y = y           Y = b      where     z = x+i.y
                                                                              X = x    to    X = a        and   f(z) = f(x+i.y) = a+i.b
                       ( Registers R03 R04 R13 R14 ... etc ... may be used by this subroutine )



001  LBL "RTZ2"
002  STO 01
003  RDN
004  STO 02
005  X<>Y
006  STO 05
007  STO 06
008  STO 09
009  STO 10
010  ST- T
011  -
012  R^
013  XEQ IND 00
014  STO 07
015  X<>Y
016  STO 08
017  RCL 02
018  RCL 01
019  RCL 05
020  ST+ X
021  ST- Z
022  -
023  XEQ IND 00
024  STO 11
025  X<>Y
026  STO 12
027  LBL 01
028  VIEW 01
029  RCL 02
030  RCL 01
031  XEQ IND 00
032  STO 03
033  X<>Y
034  STO 04
035  X<>Y
036  RCL 10
037  RCL 09
038  XEQ "Z*"
039  STO 13
040  X<>Y
041  STO 14
042  RCL 06
043  RCL 10
044  +
045  STO 20
046  RCL 05
047  RCL 09
048  +
049  STO 19
050  RCL 14
051  RCL 13
052  XEQ "Z*"
053  RCL 06
054  RCL 05
055  XEQ "Z*"
056  R-P
057  X=0?
058  GTO 02        ( a three-byte GTO )
059  CHS
060  STO 15
061  X<>Y
062  STO 16
063  RCL 06
064  RCL 20
065  +
066  RCL 05
067  RCL 19
068  +
069  RCL 14
070  RCL 13
071  XEQ "Z*"
072  STO 17
073  X<>Y
074  STO 18
075  RCL 20
076  RCL 19
077  RCL 08
078  RCL 07
079  XEQ "Z*"
080  ST- 13
081  X<>Y
082  ST- 14
083  X<>Y
084  RCL 20
085  RCL 19
086  XEQ "Z*"
087  ST- 17
088  X<>Y
089  ST- 18
090  RCL 06
091  RCL 05
092  RCL 12
093  RCL 11
094  XEQ "Z*"
095  ST+ 13
096  X<>Y
097  ST+ 14
098  X<>Y
099  RCL 06
100  RCL 05
101  XEQ "Z*"
102  ST+ 17
103  X<>Y
104  ST+ 18
105  RCL 14
106  RCL 13
107  R-P
108  STO 13
109  X<>Y
110  STO 14
111  RCL 18
112  RCL 17
113  R-P
114  RCL 16
115  ST- 14
116  ST- Z
117  X<> 15
118  ST/ 13
119  ST+ X
120  /
121  P-R
122  STO 15
123  X<>Y
124  STO 16
125  X<>Y
126  XEQ "Z^2"
127  RCL 14
128  RCL 13
129  P-R
130  XEQ "Z+"
131  XEQ "SQRTZ"
132  RCL 16
133  RCL Z
134  ST- 16
135  +
136  RCL 15
137  RCL Z
138  ST- 15
139  +
140  R-P
141  X<>Y
142  RCL 16
143  RCL 15
144  R-P
145  R^
146  X<Y?
147  RCL Y
148  X#0?
149  1/X
150  R^
151  CHS
152  X<>Y
153  P-R
154  ST+ 01
155  X<> 05
156  STO 09
157  X<>Y
158  ST+ 02
159  X<> 06
160  STO 10
161  RCL 03
162  X<> 07
163  STO 11
164  RCL 04
165  X<> 08
166  STO 12
167   E-8
168  LASTX
169  X>Y?
170  GTO 01      ( a three-byte GTO )
171  LBL 02
172  RCL 02
173  RCL 01
174  END

( 272 bytes / SIZE 021 )


      STACK        INPUTS      OUTPUTS
           Z             h             /
           Y             y0             y
           X             x0             x

( The 3 initial guesses are actually   z0 = x0+i.y0 ;   z0-(1+i).h ; z0-2.(1+i).h  )

Example:     With the same equation     sinh z + z2 + Pi = 0.     We load a routine to compute f(z)

  LBL "T"
  STO 03
  X<>Y
  STO 04
  X<>Y
  XEQ "SHZ"
  RCL 04
  RCL 03
  XEQ "Z^2"
  XEQ "Z+"
  PI
  +
  RTN

And if we choose   z0 = 1+i  and  h = 1

   T  ASTO 00
   1  ENTER^   ENTER^
       XEQ "RTZ2"  the succesive x-values are displayed and we get   -0.278189857                  ( execution time = 156 seconds )
                                                                                          X<>Y      1.812880366

Whence    z = -0.278189857 + 1.812880366.i  

Note:     In this example,  "ROOTZ"  is faster than  "RTZ2"   but if  f is a very complicated function, the advantage is clearly with "RTZ2"


        c) The Newton's Method


-Only one initial guess  z0 = x0+i.y0  is needed but we have to compute the first derivative f '(z)
-The formula is   zk+1 = zk - f(zk)/f '(zk)

Data Registers:   • R00 = f name       • R03 = f ' name    ( these 2 registers are to be initialized before executing "NWTZ" )
                                 R01 = x  
                                 R02 = y   R04 = | f '(z) |   R05: temp
Flag: /
Subroutines:  A program which computes  f(z)   assuming  x is in X-register and y is in Y-register upon entry
               and   ----------------------------  f '(z) -----------------------------------------------------------


01  LBL "NWTZ"
02  STO 01
03  X<>Y
04  STO 02
05  LBL 01
06  VIEW 01
07  RCL 02
08  RCL 01
09  XEQ IND 03
10  R-P
11  STO 04
12  X<>Y
13  STO 05
14  RCL 02
15  RCL 01
16  XEQ IND 00
17  R-P
18  ENTER^
19  X<> 04
20  X#0?
21  /
22  X<>Y
23  RCL 05
24  -
25  X<>Y
26  P-R
27  ST- 01
28  X<>Y
29  ST- 02
30   E-8
31  LASTX
32  X>Y?
33  GTO 01
34  RCL 04
35  RCL 02
36  RCL 01
37  END

( 55 bytes / SIZE 006 )


      STACK        INPUTS      OUTPUTS
           Z             /          | f(z) |
           Y             y0             y
           X             x0             x

 -The output in Z-register should be a "small" number  if  z = x+i.y  is a "true" root

Example:       sinh z + z2 + Pi = 0.     We have to load 2 routines to compute f(z) and f ' (z) = cosh z + 2.z  for instance:

  LBL "T"
  STO 07
  X<>Y
  STO 08
  X<>Y
  XEQ "SHZ"
  RCL 08
  RCL 07
  XEQ "Z^2"
  XEQ "Z+"
  PI
  +
  RTN
  LBL "U"
  STO 07
  X<>Y
  STO 08
  X<>Y
  XEQ "CHZ"
  RCL 08
  ST+ X
  RCL 07
  ST+ X
  XEQ "Z+"
  RTN

Then    T  ASTO 00   U  ASTO 03   and with  z0 = 1+i
           1  ENTER^    XEQ "NWTZ"  the successive  x-approximations are displayed
                                      and we get  -0.278189857
                                              RDN    1.812880366     whence  z = -0.278189857 + 1.812880366.i      ( execution time = 54 seconds )

Note:   | f(z) |  in Z-register actually concerns the next to last approximation.


        d) Another method

-Starting with one initial guess   z0 = x0+i.y0  , we now use the first and second derivatives   f '(z)  and  f "(z)  ,  the recursion formula is:
   zk+1 = zk -  1 / [ f '(zk)/f(zk) + E.( ( f(zk)/(2.f '(zk)) )2 - f "(zk) /(2.f (zk)) )1/2 ]   where   E = +1 or -1   ( the sign is choosen to produce the larger denominator )

Data Registers:   • R00 = f name       • R03 = f ' name     • R04 = f " name  ( these 3 registers are to be initialized before executing "RTZ4" )
                                 R01 = x  
                                 R02 = y   R05 = | f '(z) |    R06 thru R09: temp
Flag: /
Subroutines:  A program which computes  f(z)   assuming  x is in X-register and y is in Y-register upon entry
                        ----------------------------  f '(z) -----------------------------------------------------------
               and   ----------------------------  f "(z) -----------------------------------------------------------


01  LBL "RTZ4"
02  STO 01
03  X<>Y
04  STO 02
05  LBL 01
06  VIEW 01
07  RCL 02
08  RCL 01
09  XEQ IND 03
10  R-P
11  STO 06
12  X<>Y
13  STO 07
14  RCL 02
15  RCL 01
16  XEQ IND 04
17  R-P
18  STO 08
19  X<>Y
20  STO 09
21  RCL 02
22  RCL 01
23  XEQ IND 00
24  R-P
25  STO 05
26  X=0?
27  GTO 02
28  ST+ X
29  ST/ 06
30  ST/ 08
31  X<>Y
32  ST- 07
33  ST- 09
34  RCL 07
35  ST+ X
36  RCL 06
37  X^2
38  P-R
39  RCL 09
40  RCL 08
41  P-R
42  ST- Z
43  RDN
44  ST- Z
45  RDN
46  R-P
47  SQRT
48  X<>Y
49  2
50  /
51  X<>Y
52  P-R
53  RCL 07
54  RCL 06
55  P-R
56  STO 06
57  RDN
58  STO 07
59  RCL Z
60  ST- 07
61  +
62  RCL 06
63  RCL Z
64  ST- 06
65  +
66  R-P
67  X<>Y
68  RCL 07
69  RCL 06
70  R-P
71  R^
72  X<Y?
73  RCL Y
74  X#0?
75  1/X
76  R^
77  CHS
78  X<>Y
79  P-R
80  ST- 01
81  X<>Y
82  ST- 02
83   E-8
84  LASTX
85  X>Y?
86  GTO 01
87  LBL 02
88  RCL 05
89  RCL 02
90  RCL 01
91  END

( 123 bytes / SIZE 010 )



      STACK        INPUTS      OUTPUTS
           Z             /          | f(z) |
           Y             y0             y
           X             x0             x

 -The output in Z-register should be a "small" number  if  z = x+i.y  is a "true" root

Example:       sinh z + z2 + Pi = 0.     We have to load 3 routines to compute f(z) ; f ' (z) = cosh z + 2.z and  f "(z) = sinh z + 2  ,  for instance:

   LBL "T"
   STO 10
   X<>Y
   STO 11
   X<>Y
   XEQ "SHZ"
   RCL 11
   RCL 10
   XEQ "Z^2"
   XEQ "Z+"
   PI
   +
   RTN
   LBL "U"
   STO 10
   X<>Y
   STO 11
   X<>Y
   XEQ "CHZ"
   RCL 11
   ST+ X
   RCL 10
   ST+ X
   XEQ "Z+"
   RTN
   LBL "V"
   XEQ "SHZ"
   2
   +
   RTN


Then    T  ASTO 00   U  ASTO 03   V ASTO 04  and with  z0 = 1+i
           1  ENTER^    XEQ "RTZ4"    the successive  x-approximations are displayed
                                      and we get  -0.278189857
                                              RDN    1.812880366     whence  z = -0.278189857 + 1.812880366.i      ( execution time = 70 seconds )

Notes:  1-   | f(z) |  in Z-register actually concerns the next to last approximation.
             2- Among these 4 root-finding programs, "RTZ4" has the best rate of convergence.

Go back to the HP-41 software library
Go back to the general software library
Go back to the main exhibit hall