This program is Copyright © 2005 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
-Ramanujan numbers Tau(n) are defined by
x.[ (1-x).(1-x2).(1-x3)....... ]24 = Tau(1).x + Tau(2).x2 + Tau(3).x3 + ..... + Tau(n).xn + ......
-"TAU" calculates Tau(n) by the recurrence relation:
(n-1).Tau(n) = SUM1<=m<=b(n) (-1)m+1 (2m+1).(n-1-9m(m+1)/2).Tau(n-m(m+1)/2) where b(n) = INT ((8n+1)1/2-1)/2
-However, this algorithm is very sensitive to roundoff errors and results
are exact up to n = 43 only.
-Therefore, another approach is necessary to compute Tau(n) for larger
arguments.
-"TAU2" employs the formula:
Tau(n) = n4.s(n) - 24 SUM0<k<n k2.(35.k2-52.k.n+18.n2).s(k).s(n-k) ( n > 1 ) where s(k) is the sum of the divisors of k.
-Multiple precision is used and Tau(n) may be exactly computed up to
n = 13868.
1°) Program#1 ( small arguments )
-This program calculates and stores Tau(1) , Tau(2) , ...... , Tau(n) into registers R01 , R02 , ....... , Rnn
Data Registers: R00 = 0 , R01 = Tau(1) ,
........... , Rnn = Tau(n)
Flags: /
Subroutines: /
-Synthetic registers M , N , O may be replaced by any unused data registers,
for instance R97 , R98 , R99 if n < 97
01 LBL "TAU"
02 STO O
03 .1
04 %
05 2
06 +
07 STO M
08 CLX
09 STO 00
10 SIGN
11 STO 01
12 LBL 01
13 RCL M
14 INT
15 8
16 *
17 RCL 01
18 +
19 SQRT
20 RCL 01
21 -
22 2
23 /
24 INT
25 STO N
26 CLX
27 FIX 0
28 LBL 02
29 RCL M
30 RCL 01
31 RCL N
32 ST+ Y
33 *
34 2
35 /
36 -
37 X<>Y
38 RCL IND Y
39 LASTX
40 9
41 *
42 RCL 01
43 +
44 RCL M
45 INT
46 -
47 *
48 RCL N
49 ST+ X
50 RCL 01
51 +
52 *
53 RCL 01
54 CHS
55 RCL N
56 Y^X
57 *
58 +
59 DSE N
60 GTO 02
61 RCL M
62 INT
63 RCL 01
64 -
65 /
66 RND
67 STO IND M
68 ISG M
69 GTO 01
70 RCL O
71 SIGN
72 RCL IND L
73 CLA
74 FIX 4
75 END
( 106 bytes / SIZE nnn+1, but at least 003 )
STACK | INPUTS | OUTPUTS |
X | n | tau(n) |
L | / | n |
Example: Evaluate Tau(10)
10 XEQ "TAU" >>>> Tau(10) = -115920 in X-register and in R10 ( execution time = 32 seconds )
-We also have 1 , -24 , 252 , -1472
, 4830 , -6048 , -16744 , 84480 , -113643 ( i-e Tau(1) , ..........
, Tau(9) )
in registers R01 , R02 , R03 ,
R04 , R05 , R06 ,
R07 , R08 , R09
respectively
Notes: -Without lines 27 and 66 ( FIX
0 and RND ) , Tau(33) would already be wrong!
-This program yields Tau(84) = 6,188,510,373 whereas the correct value
is Tau(84) = 6,211,086,336
and 300 XEQ "TAU" produces a completely
meaningless result!
-Obviously, roundoff errors overwhelm the wanted function.
2°) Program#2 ( 1 < n < 13869
)
-The following program computes Tau(n) by groups of 4 digits in registers
R15 R16 R17 R18 R19 R20
Data Registers:
R01 = n ; R00 and R02 thru R20 are used for temporary
data storage and when the program stops:
R03 = Tau(n) = X-register ( approximate value )
and R15 to R20 contain the digits of Tau(n) ( exact value )
Flags: /
Subroutines: /
01 LBL "TAU2"
02 STO 00
03 STO 01
04 E4
05 STO 02
06 15.02
07 STO 14
08 CLRGX
if you don't have an HP-41CX , replace line 08 by 0 LBL 00
STO IND Y ISG Y GTO 00 or add CLRG
after line 01
09 DSE 00
10 LBL 01
11 RCL 00
12 X^2
13 ENTER^
14 STO Z
15 RCL 02
16 MOD
17 STO 04
18 -
19 RCL 02
20 /
21 STO 03
22 CLX
23 35
24 *
25 RCL 00
26 RCL 01
27 *
28 52
29 *
30 -
31 RCL 01
32 X^2
33 18
34 *
35 +
36 STO Y
37 RCL 02
38 MOD
39 STO 07
40 -
41 RCL 02
42 /
43 STO Y
44 RCL 02
45 MOD
46 STO 06
47 -
48 RCL 02
49 /
50 STO 05
51 RCL 00
52 XEQ 05
53 RCL 01
54 RCL 00
55 -
56 XEQ 05
57 *
58 STO Y
59 RCL 02
60 MOD
61 STO 09
62 -
63 RCL 02
64 /
65 STO 08
66 XEQ 08
67 DSE 00
68 GTO 01
69 RCL 14
70 24
71 CHS
72 LBL 02
73 ST* IND Y
74 ISG Y
75 GTO 02
76 RCL 01
77 X^2
78 STO Y
79 RCL 02
80 MOD
81 STO 04
82 STO 07
83 -
84 RCL 02
85 /
86 STO 03
87 STO 06
88 CLX
89 STO 05
90 RCL 01
91 XEQ 05
92 STO Y
93 RCL 02
94 MOD
95 STO 09
96 -
97 RCL 02
98 /
99 STO 08
100 XEQ 08
101 RCL 02
102 RCL 02
103 RCL 14
104 STO 00
105 CLX
106 LBL 03
107 RCL IND 00
108 +
109 *
110 ISG 00
111 GTO 03
112 LASTX
113 STO 03
114 SIGN
115 STO 04
116 ST* 02
117 20.015
118 STO 00
119 CLX
120 LBL 04
121 RCL IND 00
122 +
123 STO Y
124 RCL 02
125 MOD
126 STO IND 00
127 -
128 RCL 02
129 /
130 RCL 04
131 *
132 DSE 00
133 GTO 04
134 ST+ 15
135 RCL 14
136 RCL 03
137 RTN
138 LBL 05
lines 138 to 161 calculates the sum of the divisors of the number in X-register.
139 STO 08
140 SQRT
141 INT
142 STO 09
143 CLX
144 LBL 06
145 RCL 08
146 RCL 09
147 /
148 FRC
149 X#0?
150 GTO 07
151 X<> L
152 ST+ Y
153 RCL 09
154 X#Y?
155 ST+ Z
156 RDN
157 LBL 07
158 RDN
159 DSE 09
160 GTO 06
161 RTN
162 LBL 08
lines 162 to 207 calculates the product P = k2.(35.k2-52.k.n+18.n2)
and stores the result in registers R10 thru R13
163 RCL 04
164 RCL 07
165 *
166 STO Y
167 RCL 02
168 MOD
169 STO 13
170 -
171 RCL 02
172 /
173 RCL 03
174 RCL 07
175 *
176 +
177 RCL 04
178 RCL 06
179 *
180 +
181 STO Y
182 RCL 02
183 MOD
184 STO 12
185 -
186 RCL 02
187 /
188 RCL 03
189 RCL 06
190 *
191 +
192 RCL 04
193 RCL 05
194 *
195 +
196 STO Y
197 RCL 02
198 MOD
199 STO 11
200 -
201 RCL 02
202 /
203 RCL 03
204 RCL 05
205 *
206 +
207 STO 10
208 RCL 09
lines 208 to 275 computes the product of P by s(k).s(n-k) and add the result
to the partial sum in registers R15 thru R20.
209 RCL 13
210 *
211 STO Y
212 RCL 02
213 MOD
214 ST+ 20
215 -
216 RCL 02
217 /
218 RCL 08
219 RCL 13
220 *
221 +
222 RCL 09
223 RCL 12
224 *
225 +
226 STO Y
227 RCL 02
228 MOD
229 ST+ 19
230 -
231 RCL 02
232 /
233 RCL 08
234 RCL 12
235 *
236 +
237 RCL 09
238 RCL 11
239 *
240 +
241 STO Y
242 RCL 02
243 MOD
244 ST+ 18
245 -
246 RCL 02
247 /
248 RCL 08
249 RCL 11
250 *
251 +
252 RCL 09
253 RCL 10
254 *
255 +
256 STO Y
257 RCL 02
258 MOD
259 ST+ 17
260 -
261 RCL 02
262 /
263 RCL 08
264 RCL 10
265 *
266 +
267 STO Y
268 RCL 02
269 MOD
270 ST+ 16
271 -
272 RCL 02
273 /
274 ST+ 15
275 RTN
276 END
( 352 bytes / SIZE 021 )
STACK | INPUTS | OUTPUTS |
Y | / | 15.020 |
X | n | tau(n) ( approx.) |
15.020 = control number of the solution ( in registers R15 thru R20 )
Example1: Calculate Tau(84)
84 XEQ "TAU2" >>>> Tau(84) = 6211086336 ( 13m07s ) In this case, the approximate value is quite exact!
which is confirmed by R15 = 0 , R16 = 0 , R17 = 0 , R18 = 62 , R19 = 1108 , R20 = 6336
Example2: Compute Tau(314)
314 XEQ "TAU2" >>>> Tau(314) = -3.156280211 1013 = approximate value ( in 1h07mn20s )
and we have: R15 = 0 , R16 = 0 , R17 = -31 , R18 = -5628 , R19 = -210 ( which is to be read -0210 ) , R20 = -5744
whence Tau(314) = -31562802105744
-Remark that registers R15 to R20 are inferior or equal to zero when Tau(n) < 0
Notes: -With a good 41-emulator in turbo mode 13868 XEQ "TAU2" yields:
R15 = -336 ; R16 = -1948 ; R17 = -0538 ; R18 = -4247 ; R19 = -3535 ; R20 = -1552 ( in about 1h23mn on my PC )
whence Tau(13868) = -33619480538424735351552 with a "real" HP-41, execution time is probably of the order of 8 or 9 days!
-For n > 13868 , the calculations involve numbers greater than 1010
and the results are only approximate,
but if we use more registers to perform multiplications and additions,
this limit can be overcome.
Go back to the HP-41 software library
Go back to the general software library
Go
back to the main exhibit hall