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