More 41 birthday fun.



#2

Ladies, Gents, and others, here is a simple programming challenge:

What are my lottery numbers?

To find my lottery numbers you just need to find all the positive integers such that x2 - x + n is a prime number for x = 0, 1, ..., n - 1. (NOTE: 1 is not a prime number.)

Rules:

  1. 41C, 41CV, 41CX RPN only.
  2. No modules.
  3. Enter as many times as you like. I like community/group solutions where speed and size gets optimized over time.

Prizes:

Well, if I win the lottery, then the authors of the shortest and the fastest programs will get a cut. The lottery drawing is July 31st.

Edit: typo


Edited: 24 July 2009, 6:00 p.m. after one or more responses were posted


#3

Is that for n in 0,1, ..., x-1? And is that for all n in that range?


#4

Oops! I see the typo, fixed. To be clear, here is an example:

Take n = 2, therefore x = 0, 1:

x = 0, 0 - 0 + 2, prime, passed
x = 1, 1 - 1 + 2, prime, passed
All passed. 2 is one of the 6 numbers. What are the other 5?

Take n = 7, therefore x = 0, 1, ..., 6:

x = 0, 0 - 0 + 7, prime, passed
x = 1, 1 - 1 + 7, prime, passed
x = 2, 4 - 2 + 7, not prime, fails
If any one test fails, then n is not a solution.

Edited: 24 July 2009, 6:15 p.m.


#5

well; 2,3,5 and 11 are apparent without a calculator.

#6

Well, I don't have my 41 here, but I get the first 4 numbers with paper and pencil easily. Time to go to bed ;)

Edit: Ooops, Hugh was faster than me.

Edited: 24 July 2009, 6:50 p.m.


#7

now, i remember where i've seen that polynomial before. turns out it's not easy to prove there are only 6 solutions in general (although 49 is a given max for the problem here).

looks like the easiest way to prime test is trial division because the numbers are small.

#8

Well, neither the shortest nor the fastest, but it's a start.

I haven't tried it on the 41CX yet. This was done on Free42.

Gerson.

00 { 69-Byte Prgm }
01LBL "HP41BDF"
02 STO 04
03 2
04 STO 05
05LBL 01
06 RCL 05
07 1
08 -
09 STO 06
10 1
11 STO 07
12LBL 02
13 RCL 07
14 ENTER
15 1
16 -
17 *
18 RCL+ 05
19 XEQ "ISPRIM?"
20 X=0?
21 GTO 03
22 1
23 STO+ 07
24 RCL 06
25 RCL 07
26 X<=Y?
27 GTO 02
28 RCL 05
29 PRX
30LBL 03
31 1
32 STO+ 05
33 RCL 04
34 RCL 05
35 X<=Y?
36 GTO 01
37 RTN
38 .END.

00 { 70-Byte Prgm }
01LBL "ISPRIM?"
02 3
03 X<>Y
04 X>=Y?
05 GTO 00
06 1/X ; Lines 06 and 07 by Tony Duell - On the 42S they may be replaced by "3 XOR"
07 STO+ ST X ; http://www.hpmuseum.org/cgi-sys/cgiwrap/hpmuseum/archv015.cgi?read=74095 (see message #4)
08LBL 00
09 STO 00
10 SQRT
11 STO 01
12 0
13 STO 02
14 2
15 STO 03
16LBL 01
17 RCL 00
18 RCL/ 03
19 IP
20 RCL 00
21 1
22 -
23 RCL/ 03
24 IP
25 -
26 STO+ 02
27 1
28 STO+ 03
29 RCL 01
30 RCL 03
31 X<=Y?
32 GTO 01
33 RCL 02
34 +/-
35 RCL 00
36 /
37 XEQ "INT"
38 1
39 +
40 RTN
41 END

00 { 19-Byte Prgm }
01LBL "INT"
02 ENTER
03 IP
04 X<=Y?
05 GTO 01
06 1
07 -
08LBL 01
09 X<>Y
10 Rv
11 RTN
12 END

50 XEQ HP41BDF ->

2 ***
3 ***
5 ***
11 ***
17 ***
41 ***

IsPrime(x) function by Sebastián Martín-Ruiz:

http://en.wikipedia.org/wiki/Formula_for_primes


Edited: 24 July 2009, 11:43 p.m.


#9

I know this is meant to be 41 RPN only but I thought I'd submit a RPL version just in case it has any interest. It uses ISPRIME? so only works on any CAS equiped RPL machine. Note that I always work with a 48GX so don't have to deal with the exact/approx mode in the later machines. Consequently, there might be flags that need setting - putting it into approx mode for example - and I am not familiar with testing in this manner (2. <> 2) etc.. Note that the SEQ and SigmaLIST commands can probably be combined into a single Sigma command.

I accept that performance wise, this is not the most efficient as it evaluates the entire sequence for each value of n but using SEQ keeps the code simpler and this program is a one-off anyway.


\<< {} 2 \-> r n

\<<

DO

\<< x SQ x - n + ISPRIME? \>> 'x' 0 n 1 - 1 SEQ

SigmaLIST

n ==

\<< n 'r' STO+ \>> IFT

'n' INCR DROP

UNTIL

r SIZE 6 ==

END

r

\>>

\>>
The program returns a list with the lottery numbers in it as:

{ 41 17 11 5 3 2 }

If any RPL people can suggest ways of improving the code, I'd be glad to hear it!

Mark


Edited: 25 July 2009, 8:46 a.m.


#10

I whipped up an RPL version too, just for the hell of it. It's pretty verbose and unoptimized, but here we go:

\<< 0 1 0 1 \-> t n c p
\<<
WHILE t 6 <
REPEAT 0 'c' STO 1
'p' STO
WHILE p c n < AND
REPEAT c n \-> X N
'SQ(X)-X+N' ISPRIME? 'p'
STO 'c' INCR DROP
END
IF p
THEN n DUP 1 DISP
't' INCR DROP
END 'n' INCR DROP
END t \->LIST
\>>
\>>

I did it on my 48GX, so I had to make my own ISPRIME?. It's a fairly naive algorithm, but it gives the correct results for -1, 0, 1, and 2.

\<< ABS 0 1 2 \-> X SX p n
\<<
CASE X 1 ==
THEN 0
END X 2 ==
THEN 1
END X 2 MOD NOT
THEN 0
END X \SQRT IP 'SX'
STO
WHILE p n SX \<= AND
REPEAT X n MOD 'p'
STO 2 'n' STO+
END p SIGN
END
\>>
\>>

#11

  1. Your program ISPRIME? didn't work on my 48GX. I had to replace \SQRT by \v/ which translated to SQRT.

  2. When running your main program I get the following list:

    { 2 3 5 7 9 11 }

    This is not the solution I think we should get.

  3. According to your program ISPRIME? all the numbers 9, 15, 21, 27, ... are considered primes.

#12

Not that I consider myself a RPL expert but you run the code to check whether the numbers are prime over the whole sequence. But as soon as one item fails you could stop and reject that number. I'm afraid you can't do that with SEQ. An optimization would be to use only primes as n. Another would be to use only prime twins. But then you're already on your way to roll up your program.

Besides being a little slow a very elegant solution.

#13

A couple of obvious optimizations but still rather slow: 14min 09s on the HP-41CV. Definitely LOG and 10^x in the floor function routine is not a good idea, even thought they don't appear to take too long for those arguments. Also, a faster primality test would help...

Anyway, compared to the other 41 birthday challenge, this is a piece of cake :-)

01LBL'41BD
02 0
03 STO 08
04 1
05 STO 04
06 2
07 STO 05
08LBL 01
09 RCL 05
10 XEQ'ISPR?
11 X=0?
12 GTO 03
13 RCL 05
14 1
15 -
16 STO 06
17 1
18 STO 07
19LBL 02
20 RCL 07
21 ENTER^
22 1
23 -
24 *
25 RCL 05
26 +
27 XEQ'ISPR?
28 X=0?
29 GTO 03
30 1
31 ST+ 07
32 RCL 06
33 RCL 07
34 X<=Y?
35 GTO 02
36 RCL 05
37 VIEW X
38 1
39 ST+ 08
40LBL 03
41 RCL 04
42 ST+ 05
43 2
44 STO 04
45 5
46 RCL 08
47 X<=Y?
48 GTO 01
49 RTN
50 END

01LBL'ISPR?
02 3
03 X<>Y
04 X>Y?
05 GTO 00
06 1/X
07 ST+ X
08LBL 00
09 STO 00
10 SQRT
11 STO 01
12 0
13 STO 02
14 2
15 STO 03
16LBL 01
17 RCL 00
18 RCL 03
19 /
20 INT
21 RCL 00
22 1
23 -
24 RCL 03
25 /
26 INT
27 -
28 ST+ 02
29 1
30 ST+ 03
31 RCL 01
32 RCL 03
33 X<=Y?
34 GTO 01
35 RCL 02
36 CHS
37 RCL 00
38 /
39 0
40 X>Y?
41 10^X
42 10^X
43 LOG
44 -
45 INT
46 1
47 +
48 RTN
49 END

Edit to fix some typos


Edited: 25 July 2009, 9:41 p.m.


#14

Slightly less slow and shorter:

01LBL'41BD
02 2
03 VIEW X
04 1
05 STO 04
06 3
07 STO 05
08LBL 01
09 RCL 05
10 XEQ'ISPR?
11 X=0?
12 GTO 03
13 2
14 STO 06
15-LBL 02
16 RCL 06
17 X^2
18 LASTx
19 
20 RCL 05
21 +
22 XEQ'ISPR?
23 X=0?
24 GTO 03
25 1
26 ST+ 06
27 RCL 06
28 RCL 05
29 X#Y?
30 GTO 02
31 RCL 05
32 VIEW X
33 1
34 ST+ 04
35-LBL 03
36 2
37 ST+ 05
38 5
39 RCL 04
40 X<=Y?
41 GTO 01
42 END

01LBL'ISPR?
02 3
03 X<>Y
04 X>Y?
05 GTO 00
06 1/X
07 ST+ X
08LBL 00
09 STO 00
10 SQRT
11 STO 01
12 0
13 STO 02
14 2
15 STO 03
16LBL 01
17 RCL 00
18 RCL 03
19 /
20 INT
21 RCL 00
22 1
23 -
24 RCL 03
25 /
26 INT
27 -
28 ST+ 02
29 1
30 ST+ 03
31 RCL 01
32 RCL 03
33 X<=Y?
34 GTO 01
35 RCL 02
36 RCL 00
37 /
38 X<=0?
39 GTO 02
40 1
41 +
42-LBL 02
43 CHS
44 INT
45 1
46 +
47 RTN
48 END

2 00min 00s
3 00min 04s
5 00min 14s
11 01min 02s
17 02min 46s
41 12min 51s


#15


Hi, I'm glad to see an 41 implementation. I was going to point out that no one had entered an _actual_ 41 program, although i am also guilty :-)

Looking into the prime test, the smallness of the numbers means we can get away with trivial division. here's a short LUA program where i've rolled all the stages together into 3 loops,

LUA syntax is pretty easy to see what's going on. the only non-obvious syntax is n\k means floor(n/k).

function findlots2()
-- condensed into one function
for i = 2,49 do
ok = true
for j = 0, i - 1 do
n = j*j - j + i
for k = 2,math.sqrt(n) do
if n\k * k == n then ok = false; break end
end
if not ok then break end
end
if ok then print(i) end
end
end

i would think this would code up to a fairly short 41 program, especially since there is an extra test to escape the middle loop which could actually just be a goto.


#16

Your concise algorithm fits even the HP-12C! I tried it on the newest 12C+. It takes 25 seconds to run or about 25 minutes on a normal 12C (and under 10 minutes on the 41CX, I estimate). Although the Pascal code works, on the 12C 2 doesn't show up. I have to find out what's gone wrong with that.

Regards,

Gerson.

01  2                   ; for i:=2 to 49 do 
02 STO 0 ; begin
03- 0
04 STO 4 ; ok:=true
05 1 ; for j:=1 to i-1 do
06 STO 1 ; begin
07- RCL 1
08 ENTER
09 ENTER
10 1
11 -
12 *
13 RCL 0
14 +
15 STO 3 ; n:=j*(j-1)+i
16 2 ; for k:=2 to Sqrt(n) do
17 STO 2 ; begin
18- RCL 3
19 RCL 2
20 /
21 INTG
22 RCL 2
23 *
24 RCL 3
25 -
26 x=0 ; if (n div k)*k = n then
27 GTO 29 ; ok:=false
28 GTO 32
29- 1
30 STO 4
31 GTO 39 ; exit loop
32- 1
33 STO+ 2
34 RCL 3
35 SQRT
36 RCL 2
37 X<=Y ; end
38 GTO 18
39- RCL 4 ; if not ok then
40 x=0
41 GTO 43
42 GTO 51 ; exit loop
43- 1
44 STO+ 1
45 RCL 0
46 1
47 -
48 RCL 1
49 x<=y
50 GTO 07 ; end
51- RCL 4
52 X=0 ; if ok then
53 GTO 55 ;
54 GTO 57
55- RCL 0
56 PSE ; Write(i)
57- 1
58 STO+ 0
59 4
60 9
61 RCL 0
62 x<=y
63 GTO 03 ; end

Program HP41BD;
var ok: boolean;
i, j, k, n: integer;
begin
for i:=2 to 49 do
begin
ok:=true;
for j:=1 to i-1 do
begin
n:=j*j-j+i;
for k:=2 to Trunc(Sqrt(n)) do
if n div k * k = n then ok:=false;
if not ok then ;
end;
if ok then WriteLn(i:3);
end
end.

Running
2
3
5
11
17
41


#17

I like it!

congrats to you. I think you deserve to win the competition as no one else (including myself) has an actual 41 program.

it does show just how slow the real 41 is though. 25 seconds vs 10 mins. hmmm.


#18

As I was guilty of posting a non-41 RPN program, I have to confess I had an ulterior motive behind posting that RPL program.

With so many people not being fans of RPL, I wanted to show that RPL could solve this problem in a fairly elegant manner with quite a simple program and I also wanted to use lists in the solution as that is a key part of the language.

Mark


#19

I think there is no problem posting non-41 solutions. There's always a chance the algorithms may be interesting and useful. I thought of posting an RPL version solution too but I don't think it would be any better than the ones you have already posted.

Gerson.

#20

Exactly 10min 22s on the HP-41CX, but the credit goes to you because of the algorithm. Anyway, the odds of Egan's numbers winning the lottery are very dim :-)

Since 2 is very obvious, it is displayed at once, then only the odd numbers are tested (but the increase in speed should not be significant because of this). I don't know the 41 enough to try further speed optimizations though.

         2   00min 00s
3 00min 04s
5 00min 11s
11 00min 48s
17 02min 12s
41 10min 22s
I cannot post the program: I've been receiving the message "Sorry, data too long for this forum" but the program is only 68 stepts long. I'll try later.
----

01LBL "HP41BD"
02 3.04902
03 STO 00
04 2
05 VIEW X
06LBL 01
07 0
08 STO 04
09 1
10 STO 01
11LBL 02
12 RCL 01
13 ENTER
14 ENTER
15 1
16 -
17 *
18 RCL 00
19 INT
20 +
21 STO 03
22 2
23 STO 02
24LBL 03
25 RCL 03
26 RCL 02
27 /
28 INT
29 RCL 02
30 *
31 RCL 03
32 X#Y?
33 GTO 05
34 1
35 STO 04
36 GTO 06
37LBL 05
38 1
39 ST+ 02
40 RCL 03
41 SQRT
42 RCL 02
43 X<=Y?
44 GTO 03
45LBL 06
46 RCL 04
47 X#0?
48 GTO 08
49 1
50 ST+ 01
51 RCL 00
52 INT
53 1
54 -
55 RCL 01
56 X<=Y?
57 GTO 02
58LBL 08
59 RCL 04
60 X#0?
61 GTO 09
62 RCL 00
63 INT
64 VIEW X
65LBL 09
66 ISG 00
67 GTO 01
68 .END.

Edited to include program.


Edited: 27 July 2009, 10:48 p.m. after one or more responses were posted


#21

Another optimisation would be to only test primes instead of odd numbers. Don't know if this can be incorporated easily or not.

In this case this would be equivalent to testing for divisibility by 2, 3 & 5. Add 7 if your lotto includes 49. Doing the even optimisation which is already there, one 3 & 5 need to be checked.


- Pauli


#22

What really takes time are the primality tests required for x2 - x + n, for x = 1, 2, ..., n - 1. For instance, when n = 41, 40 tests have to be performed. Skipping the few non-primes in the range 2..49 would have a very low effect in the reduction of the overall running time.

Gerson.

#23

Gerson,

beautiful as always!Having done my own little code, I now have the pleasure to see the gems posted already.

Really nifty prime checking algorithm, not the clunky thing I have!

I wonder how much your speed would increase if you replace the numbers with registers. Digits really are very slow on the hp41. Another 'expensive' calculation on the 41 is SQRT so I tried to not have it in a loop. I can't quite tell if you recalc it more than once, but probably not :-)

Also, I saw that Paul suggested to only test for prime numbers as I do in my code. However I do think that makes a difference in run-time. For example, 49 gets very quickly abandoned without ever checking a single x.

Lastly, I think that starting the x2-x+n test with n-1 is advantagious. It is the largest number to be tested for prime and hence has the highest chance to be not prime, at which point no further tests have to be done. This is why I started from 49 -> 2 rather then the other way around.

Cheers

Peter

#24

Gerson,

beautiful as always! Having done my own little code, I now have the pleasure to see the gems posted already.

Really nifty prime checking algorithm, not the clunky thing I have!

I wonder how much your speed would increase if you replace the numbers with registers. Digits really are very slow on the hp41. Another 'expensive' calculation on the 41 is SQRT so I tried to not have it in a loop. I can't quite tell if you recalc it more than once, but probably not :-)

Also, I saw that Paul suggested to only test for prime numbers as I do in my code. However I do think that makes a difference in run-time. For example, 49 gets very quickly abandoned without ever checking a single x.

Lastly, I think that starting the x2-x+n test with n-1 is advantagious. It is the largest number to be tested for prime and hence has the highest chance to be not prime, at which point no further tests have to be done. This is why I started from 49 -> 2 rather then the other way around.

Cheers

Peter


#25

Quote:
Really nifty prime checking algorithm,

Thats Hugh Steer's algorithm, to whom the compliments go. I just ported his Lua code to the 41, not in an efficient way though as I've never been an HP-41 user (back in the day I had a 15C).
Your optimization tips are very much appreciated, thanks!

Gerson.

#26

Quote:
Also, I saw that Paul suggested to only test for prime numbers as I do in my code. However I do think that makes a difference in run-time. For example, 49 gets very quickly abandoned without ever checking a single x.

I did this in my first 41 program (Edited: 25 July 2009, 9:41 p.m.). I remember the running time dropped from 14m50s to 14m09s, not a significant gain to justify the extra steps, I think.

Regards,

Gerson.


Edited: 29 July 2009, 5:29 p.m.


#27

I fiddled a bit and an subroutine to do a full check for being prime doesn't seem to help any (actually, it hurts).

Just checking for divisibility by 3 and 5 ought to help a little, however the first evaluation of the polynomial will catch these out anyway so the saving will likely be slight.

Storing the first dozen or so primes in registers and using these for both the outer loop and the primality tester might be a win.

The other thing that looks bad in the inner loop is the square root which should also be avoidable. I'm thinking along the lines of

if sn = sqrt(n2 - n + p) then, sn+1 = 1 + sn
and start s0 as sqrt(p).

- Pauli

#28

I stand corrected. doing a similar test no my machine when starting with x=1 it takes longer to have the extra prime test.

#29

I like the new compact algorithm.

The improvement I would make is to only use odd numbers for the trial division. Since we're only testing odd n, and x*(x-1) + n must be odd, we don't need to divide by 2 (or other even numbers).

I only modified Gerson's code slightly, but it executed in just over 6 minutes on my HP-41 CX.

         2   00min 00s
3 00min 04s
5 00min 10s
11 00min 36s
17 01min 26s
41 06min 06s

01 LBL "HP41BD"
02 3.04902
03 STO 00
04 2
05 VIEW X
06 LBL 01
07 0
08 STO 04
09 1
10 STO 01
11 LBL 02
12 RCL 01
13 ENTER
14 ENTER
15 1
16 -
17 *
18 RCL 00
19 INT
20 +
21 STO 03
22 1
23 STO 02
24 GTO 05
25 LBL 03
26 RCL 03
27 RCL 02
28 /
29 INT
30 RCL 02
31 *
32 RCL 03
33 X#Y?
34 GTO 05
35 1
36 STO 04
37 GTO 06
38 LBL 05
39 2
40 ST+ 02
41 RCL 03
42 SQRT
43 RCL 02
44 X<=Y?
45 GTO 03
46 LBL 06
47 RCL 04
48 X#0?
49 GTO 08
50 1
51 ST+ 01
52 RCL 00
53 INT
54 1
55 -
56 RCL 01
57 X<=Y?
58 GTO 02
59 LBL 08
60 RCL 04
61 X#0?
62 GTO 09
63 RCL 00
64 INT
65 VIEW X
66 LBL 09
67 ISG 00
68 GTO 01
69 BEEP
70 .END.

(edited to correct listing to be what was in calculator)

Edited: 28 July 2009, 4:01 p.m. after one or more responses were posted


#30

Hi Steve,

Quote:
I like the new compact algorithm.

Courtesy of Hugh Steers :-)

I do appreciate you efforts in reducing the running time. The first program I submitted would take about 15 minutes to run (provided the compability issues are resolved).

Avoiding repeatedly used numbers in program lines as Peter has suggest will help and should be easy to implement.

Gerson.

Edited: 28 July 2009, 2:01 p.m.

#31

Here are my improvements.

1. Calculate the sqrt once and store for testing.

2. Get rid of the OK flag. Not needed with GTOs.

3. Clean up the logic and use numbers already on the stack.

I kept the test from low to high because as is it gives us a free prime test for n when x=1. The code is much shorter, but not as efficient as calling a prime testing subroutine twice as PeterP has shown us.

Final execution time is now 4m 26s.

Suggestions and improvements are welcome. It wouldn't exist without the ideas of Gerson and PeterP, so my thanks to them.

1  LBL "HP41BD"
2 3.04902
3 STO 00
4 2
5 VIEW X
6 LBL 01
7 1
8 STO 01
9 LBL 02
10 RCL 01
11 ENTER
12 ENTER
13 1
14 -
15 *
16 RCL 00
17 INT
18 +
19 STO 03
20 SQRT
21 STO 05
22 3
23 STO 02
24 LBL 03
25 RCL 05
26 RCL 02
27 X>Y?
28 GTO 05
29 RCL 03
30 X<>Y
31 /
32 FRC
33 X=0?
34 GTO 09
35 2
36 ST+ 02
37 GTO 03
38 LBL 05
39 1
40 ST+ 01
41 RCL 00
42 X<>Y
43 -
44 RCL 01
45 X<=Y?
46 GTO 02
47 VIEW X
48 LBL 09
49 ISG 00
50 GTO 01
51 BEEP
52 END

#32

Quote:
Final execution time is now 4m 26s.

By replacing lines 11 through 18 with the ones below you can save 6 seconds and 2 steps:

11 X^2
12 LASTx
13 -
14 RCL 00
15 INT
16 +

#33

Nice.

My gut feeling was stack operations would be faster than X^2, so I used your original method.

Empirical evidence wins! Shorter and faster is better.

Thanks for the improvement.


#34

Quote:
My gut feeling was stack operations would be faster

So was mine. The following modification saves one more step and two registers (R02 and R05) and is even faster: 3 minutes and 57 seconds to run.

01 LBL 'BD
02 3.04902
03 STO 00
04 2
05 VIEW X
06 LBL 01
07 1
08 STO 01
09 LBL 02
10 RCL 01
11 X^2
12 LASTx
13 -
14 RCL 00
15 INT
16 +
17 STO 03
18 SQRT
19 3
20 LBL 03
21 X>Y?
22 GTO 05
23 RCL 03
24 X<>Y
25 /
26 LASTx
27 X<>Y
28 FRC
29 X=0?
30 GTO 09
31 RDN
32 2
33 +
34 GTO 03
35 LBL 05
36 1
37 ST+ 01
38 RCL 00
39 X<>Y
40 -
41 RCL 01
42 X<=Y?
43 GTO 02
44 VIEW X
45 LBL 09
46 ISG 00
47 GTO 01
48 BEEP
49 END

-----

P.S.: If line 36 is changed to

36 SIGN
as per Peter's advice, the total running time drops to 3m 54s.


Edited: 29 July 2009, 3:20 p.m.


#35

Nice use of the stack instead of registers. Great idea.

I incorporated your ideas with my latest changes. Most significant was I realized the HP-41 has a MOD command. That replaces the divide and FRC sequence and saves a bunch in the inner loop.
I also got rid of a redundant RCL after LBL 02, and used register 04 to hold the constant 2 following PeterP's idea.

Since I couldn't get PeterP's code to run on my simulator, I think this is the current leader.

The code now executes in 3min 0sec.

1  LBL "41BD"
2 3.04902
3 STO 00
4 2
5 VIEW X
6 STO 04
7 LBL 01
8 1
9 STO 01
10 LBL 02
11 X^2
12 LASTx
13 -
14 RCL 00
15 INT
16 +
17 STO 03
18 SQRT
19 3
20 LBL 03
21 X>Y?
22 GTO 05
23 RCL 03
24 X<>Y
25 MOD
26 X=0?
27 GTO 09
28 RDN
29 LASTx
30 RCL 04
31 +
32 GTO 03
33 LBL 05
34 SIGN
35 ST+ 01
36 RCL 00
37 X<>Y
38 -
39 RCL 01
40 X<=Y?
41 GTO 02
42 VIEW X
43 LBL 09
44 ISG 00
45 GTO 01
46 BEEP
47 END

#36

Great improvements, congratulations!

It appears we can break the 3-minute barrier if we change line 28 to:

28 X<>Y

2min 55sec

I've been using a real 41CX (SN 2819S21701). Notice both parts of the serial number are primes! Even S is prime, being the 19th letter :-)

Gerson.


Edited: 29 July 2009, 6:39 p.m.


#37

I was able to shave a second or two by replacing lines 28 and 29 in the inner loop with

28 X<> L

The code is also now only 46 lines long. I've also saved a little bit with a change to the code at LBL 05, but it makes little overall difference. LBL 02 has to move before the STO 01:

8  1
9 LBL 02
10 STO 01
...
32 LBL 05
33 SIGN
34 RCL 01
35 +
36 RCL 00
37 INT
38 X<>Y
39 X#Y?
40 GTO 02
41 VIEW X
42 LBL 09
...

I hope this makes sense. My calculator may be a bit slower than yours. The SN is 2850S22061, so it's clearly not optimized for prime testing.

We could also save a little bit using synthetic registers, but I think it's nicer to use instructions that any HP-41 user can enter and try. I've relearned a lot about this nice machine in the last few days.


#38

I love this collaboration! Using all the great improvements posted and incorporated I made one more change to the above code. Namely starting from x = n-1 rather than x = 1. I also added a couple of lines to replace the number 3 with a Rcl. Last but not least I used synthetic registers.

Using Steve's run-times as a benchmark I would think that the code below should execute on his machine in 61 seconds. (The last number, 41, should appear after about 51 seconds as it takes another 10 to check that there is no number between 41 and 49)

I wonder if we can drop below the 1min barrier?

Cheers

Peter

LBL EFv3
3.04902
STO 00 ;(n)
INT
STO a ;(3)
2
VIEW X
STO M ;(2)
STO N ;(x)
LBL 02 ;Loop over x
RCL N ;(x)
x^2
LastX
-
RCL 00 ;(n)
INT
+
STO O ;(x^2-x+n)
SQRT
RCL a ;(3)
LBL 03 ;primeality test
X>Y?
GTO 05 ;->it is prime
RCL O ;(x^2-x+n)
X<>Y
MOD
x=0?
GTO 09 ;->not prime
X<>L
RCL M ;(2)
+
GTO 03
LBL 05 ;next x
DSE N ;(x)
GTO 02
RCL 00 ;(n)
INT
View X ;all were prime, show result
STO N ;(prep next max x)
SIGN
ST- N
LBL 09 ;loop over n
ISG 00
GTO 02
BEEP
END

46 lines

Edited: 30 July 2009, 3:16 p.m.


#39

Wow!

That's impressive. I clock the code at 1m 04s on my emulator (which set at lowest speed matches my HP-41CX as close as my watch and fingers can tell).

I'm using the wonderful V41 emulator by Warren Furlow for these sythetic ones. I use the HP41UC User Code file converter by Leo Duran to compile the text to RAW format.

I check my own code on my real HP-41CX, but I can't do synthetics on it right now. I need to get my card reader working with a good rechargable battery for that.

I'm still looking at what you did to try to understand it, but I'm impressed.

#40

I think there is a bug in your code. Can code be buggy that produces the proper output? Hmmm.

I guess we could write code that merely displays the numbers we know we want with no calculations, but that wouldn't satisfy the requirement, right?

The problem seems to be that register N only gets set to n-1 when we successfully find and show the values that pass. It should also be set to the new n-1 when we find a non-prime and increment n for the next set of trials.

I made a fix and your code is still the fastest we have seen so far. 2m 41s.

I've included my changes to the end of your code. I really like the way you indent and comment to make the code readable, by the way. Did I missunderstand, or do you have a different idea how this should go?

 LBL 05 ;next x
DSE N ;(x)
GTO 02
RCL 00 ;(n)
INT
View X ;all were prime, show result
LBL 09 ;loop over n
RCL 00 ;(n)
INT
STO N ;(prep next max x)
ISG 00
GTO 02
BEEP
END

#41

yes, indeed you are right, my bad! Great fix!

Just for kicks and giggles, we can shave off 1 sec by replacing 'BEEP' with a synthetic tone 7 (BEEP takes a bit over a second and tone 7 about 1/10th of a second)

I wonder if there is a way to get rid of the repeated RCL 00, INT sequence. INT is not very expensive, but somewhat... But then we have to find an elegant way to do the ISG...

Cheers

Peter

Edited: 30 July 2009, 4:39 p.m.


#42

Use a second register to mirror the integer part of the counter?

Probably more fuss that it is worth though.


- Pauli

#43

Now my turn to admit to a bug.

I intialized the next x for testing to n-2. It should start at n-1.

This equation is quite forgiving. For n which don't pass, it usually fails for lots of values of x!

#44

Quote:
I was able to shave a second or two

Actually 4 and a half: now it runs in 2m 50.5s. Yet another second is gained when changing line 08 to

08 SIGN

Now, changing the outer loop, as Peter has suggested:

01 LBL 'HP41BD2
02 2
03 VIEW X
04 STO 04
05 49
06 STO 02
07 3
08 STO 00
09 LBL 01
10 SIGN
11 LBL 02
12 STO 01
13 X^2
14 LASTx
15 -
16 RCL 00
17 +
18 STO 03
19 SQRT
20 3
21 LBL 03
22 X>Y?
23 GTO 05
24 RCL 03
25 X<>Y
26 MOD
27 X=0?
28 GTO 09
29 X<>L
30 RCL 04
31 +
32 GTO 03
33 LBL 05
34 SIGN
35 RCL 01
36 +
37 RCL 00
38 X<>Y
39 X#Y?
40 GTO 02
41 VIEW X
42 LBL 09
43 RCL 04
44 ST+ 00
45 RCL 00
46 RCL 02
47 X>Y?
48 GTO 01
49 BEEP
50 END

-> 2m 47s

Or, if we want to save one step:

45 RCL 00
46 X<=NN?
47 GTO 01
48 BEEP
49 END

But then the running time rises to 2m 49s...

Just smalls steps though. Looking forward to the giant leaps you, Peter and others have been making :-)


#45

Very nice routine.

I was working on removing the ISG myself, in order to avoid all (2 here) of the INT instructions, especially in the inner loop. I couldn't get it to work as well as I hoped in the backwards testing case following Peter's code.

There is one quibble I'd make with your implementation. It looks like you're skipping the test for n=49, since in step 47 you use the X>Y test with 49 in X so all odd values stictly less than 49 get tested. This is different than what we have been doing with the ISG code, so I don't think it was intentional.

In my testing using 51 in register 02 instead of 49 only added 1 second to the code. You could also change the test to X#Y? with the 51 value, but this is the outer loop, so it only saves 24 * 14ms or 1/3 second. Yes, I did obtain a table of HP-41 function times. I'm probably working way too hard for small incremental gains.

To save the last little bit of time we could save the constant 3 in register 05, and RCL 05 at your line 20 to save a little more time but make the code 51 steps.

I'm now starting to be confused about which registers are used for which constants. Just to confound the issue, I've changed so the following apply:

REG  Usage
00 trial n
01 trial x
02 2
03 3
04 p = x^2 - x + n
05 51

Acording to my HP-41CX and my fingers on my wristwatch it takes 2m 44s. I guess I need to add stopwatch code to be more precise.

Here's my code as it currently stands:

1  LBL "41BD2"
2 2
3 VIEW X
4 STO 02
5 51
6 STO 05
7 3
8 STO 00
9 STO 03
10 LBL 01
11 SIGN
12 LBL 02
13 STO 01
14 X^2
15 LASTX
16 -
17 RCL 00
18 +
19 STO 04
20 SQRT
21 RCL 03
22 LBL 03
23 X>Y?
24 GTO 05
25 RCL 04
26 X<>Y
27 MOD
28 X=0?
29 GTO 09
30 X<> L
31 RCL 02
32 +
33 GTO 03
34 LBL 05
35 SIGN
36 RCL 01
37 +
38 RCL 00
39 X<>Y
40 X#Y?
41 GTO 02
42 VIEW X
43 LBL 09
44 RCL 02
45 ST+ 00
46 RCL 00
47 RCL 05
48 X#Y?
49 GTO 01
50 BEEP
51 END

P.S. Using the internal Stopwatch, I get 2m 44.5s (stopping before the beep).


Edited: 31 July 2009, 12:23 p.m.


#46

Quote:
There is one quibble I'd make with your implementation. It looks like you're skipping the test for n=49, since in step 47 you use the X>Y test with 49 in X so all odd values stictly less than 49 get tested. This is different than what we have been doing with the ISG code, so I don't think it was intentional.

You're right! I used to inspect register 00 after program execution to see if 51 was there, but I may have forgotten to do it after this change. My feeling was X<Y was faster than X<=Y as the test is evidently simpler.

Quote:
Using the internal Stopwatch, I get 2m 44.5s (stopping before the beep).

2m 42.5s on mine! That's 5.5 seconds faster than the previous (corrected) program. Quite another leap! Actual measurements: 2m 42.48s, 2m 42.47s, 2m 48s. Leave or take some milliseconds due to distance from ear to HP-41 buzzer. Temperature (51.8 °F) and air humidity (87 %) might affect timing also :-)

Edited: 31 July 2009, 1:24 p.m.

#47

Now it takes 2m 39.5s (2m 39.39s, 2m 39.45s; 50°F, 93% :-)
It's like trying to milk a stone...

01 LBL '41BDN
02 2
03 VIEW X
04 STO 04
05 51
06 STO 02
07 3
08 STO 00
09 STO 03
10 LBL 01
11 SIGN
12 LBL 02
13 STO 01
14 X^2
15 LASTX
16 -
17 RCL 00
18 +
19 ENTER
20 SQRT
21 RCL 03
22 LBL 03
23 X>Y?
24 GTO 05
25 RCL Z
26 X<>Y
27 MOD
28 X=0?
29 GTO 09
30 X<> L
31 RCL 04
32 +
33 GTO 03
34 LBL 05
35 SIGN
36 RCL 01
37 +
38 RCL 00
39 X<>Y
40 X#Y?
41 GTO 02
42 VIEW X
43 LBL 09
44 RCL 04
45 ST+ 00
46 RCL 00
47 RCL 02
48 X#Y?
49 GTO 01
50 BEEP
51 END


REG Usage
00 trial n
01 trial x
02 51
03 3
04 2

I know the register usage is a bit confusing but I wanted to try X#NN? (this ended up being half a second slower). Better revert to the original usage in your next improvement.


#48

I'm unfamiliar with the X#NN? instruction. Perhaps it's from an extension module that I don't have, or something? Likely something they added on the CX that I'm not familiar with.

Nice usage of the stack to save a memory register and tighten the inner loop!

I actually tried doing it backwards as Peter originally did, with the starting value given to the program in X. In my emulator it was faster, but on the real HP-41CX it wasn't. It spends the vast majority of the time verifying that 41 passes. All the other values fail quickly or pass quicker (being smaller). It took nearly 2 minutes just to display 41! It saved a memory register and made the code a bit smaller, but not as fun to watch, unless you like the flying "goose".


#49

Quote:
I'm unfamiliar with the X#NN? instruction. Perhaps it's from an extension module that I don't have, or something? Likely something they added on the CX that I'm not familiar with.

From HP-41CX Quick Reference Guide, page 25:

X=NN?      \
X#NN? | Conditional: Uses contents of
X<NN? | Rnn (NN specified in Y-register)
X<=NN? | for comparison. If not true, skips
X>NN? | next program line.
X>=NN? /

#50

Thanks for the insight.

I guess the lottery is today. Perhaps we share the winnings with Peter, OK?


#51

Let's not forget Hugh, who gave us the first version of the algorithm... provided of, course, anyone doesn't come up with a faster program. But since all Egan's numbers are odd, the odds of winning are very unlikely! Anyway, I'm pleased with what I have learned here lately :-)

#52

I decided to look at your clever stack usage one last time. If we put P=x^2-x+n both in Z and T to enter the loop at LBL 03, we can use R^ instead of RCL Z. It saves a whole 8ms, but requires that we do a double ENTER to start the loop. The auto dup of register T when doing the MOD and + instructions keeps the stack the same for the next loop.

I think I'm getting too obsessed with this, but it should save about 2 seconds. Who knew a stone could be milked?

This time I included my stopwatch code. Leave it out for a non-CX.

1  LBL "41BD2"
2 CLX
3 SETSW
4 RUNSW
5 2
6 VIEW X
7 STO 02
8 51
9 STO 05
10 3
11 STO 00
12 STO 03
13 LBL 01
14 SIGN
15 LBL 02
16 STO 01
17 X^2
18 LASTX
19 -
20 RCL 00
21 +
22 ENTER
23 ENTER
24 SQRT
25 RCL 03
26 LBL 03
27 X>Y?
28 GTO 05
29 R^
30 X<>Y
31 MOD
32 X=0?
33 GTO 09
34 X<> L
35 RCL 02
36 +
37 GTO 03
38 LBL 05
39 SIGN
40 RCL 01
41 +
42 RCL 00
43 X<>Y
44 X#Y?
45 GTO 02
46 VIEW X
47 LBL 09
48 RCL 02
49 ST+ 00
50 RCL 00
51 RCL 05
52 X#Y?
53 GTO 01
54 STOPSW
55 BEEP
56 RCLSW
57 END

#53

Quote:
I think I'm getting too obsessed with this, but it should save about 2 seconds.

Indeed! 2m 37.77s according to the HP-41CX stopwatch.

Congratulations for yet another faster version!




Edited: 31 July 2009, 9:06 p.m.

#54

brilliant idea the MOD!!!

#55

Can anyone prove that there are no such numbers > 41?
If not, can a counter example be found?

I rather doubt the second is likely given the massive increase in the number of numbers to validate as things get large. However, the first would be interesting.

- Pauli

#56

late to the party I just saw this a few hours ago -> productivity went rapidly to zero to produce the code below. It is not the shortest but it seems to find the numbers in about 3' flat on a normal CX (i41cx on default speed). I must have made some mistake...

Here are my thoughts about the code:

  • as 12 -1 +n has to be a prime, n has to be a prime. This is used by only checking odd n and doing a primality test on n itself first. If that test fails, that n is discarded and the next n = n-2 is checked
  • as n is odd, x2-x + n is also always odd. So the primality test does not check for 2
  • The primality test uses division by 3 and 5 directly and then starts by checking divisors <= sqrt(number_to_be_checked). Speed could improve if one were to start with smaller divisors towards sqrt(n), but I have not implemented this.
  • as digits are the slowest things on the 41, often used numbers are stored in registers. To give you a sense how incredibly slow numbers are, the sequence clx, sign is faster than '1'! And times for numbers are linear with the digits
  • Status registers (incl M,N,O) are used as they are slightly faster than normal registers.
  • last but not least, 2 is simply displayed as the last result, given that it was already provided as a solution

Usage: enter the highest possible number, 49, and press R/S. In about 3 minutes 02 seconds the sequence 41, 17, 11, 5, 3, 2 is displayed.

Registers:

  • STO N: number n currently tested
  • STO M: x currently tested, x e[1,n-1]
  • STO O: 1
  • STO 00: 2
  • STO 01: 3
  • STO 02: 5
  • STO 03: current divisor
  • STO 04: number tested for primality (either n or x2-x+n)

Sub-routines

  • Lbl 00: Main loop fro testing n, starting with the value the user entered to start the program
  • LBL 01: Loop to check for all x that x2-x+n is prime
  • Lbl 05: primality test
  • Lbl 06: show n as a valid result
  • Lbl 02: move to next n -> current n-2

Code:

LBL 'EF3'
STO N ;n
SIGN
STO O ;1
Enter
+
Sto 00 ;2
LastX
+
Sto 01 ;3
Rcl 00
+
Sto 02 ;5
Lbl 00
Rcl N ;n
Xeq 05 ;primality test
Rcl N ;we only return here if prime
Rcl O ;1
-
Sto M
Lbl 01 ;loop to check x^2-x-n
Rcl M
X^2
LastX
-
Rcl N
+
Xeq 05
DSE M
Gto 01 ;check next x
Gto 06 ;all prime, show n

Lbl 05 ;primality test
Sto 04
Rcl 01 ;3 no need to check 2 as only odd numbers here
X=y?
Rtn ;is prime
/
FRC
X=0? ;is not prime?
GTO 02 ;next prime
Rcl 04
Rcl 02 ;5
X=y?
RTN
/
FRC
X=0?
GTO 02
RCL 04
SQRT ;find max divisor
INT ;make sure it is odd
Rcl 00
/
INT
Rcl 00
*
RCL 0
+
Sto 03 ;divisor
Lbl 03
Rcl 02 ;make sure it is > 5
X>=Y?
RTN ; is prime
RCL 04
Rcl Z ;divisor
/
FRC
X=0?
GTO 02; not prime
Rcl 03
Rcl 00 ;2
-
Sto 03 ;next divisor
Gto 03

Lbl 06 ;show n
View N
Lbl 02 ;next n to test
Rcl N
Rcl O
-
Sto N
Rcl 01 ;is it >=3?
X<=Y?
Gto 00 ;yes, test next n
View 00 ;no. show 2 and stop
Stop

Edited: 28 July 2009, 2:13 p.m.


#57

A very nice entry, and I like your ideas.

Especially sweet as you apparently didn't look at what others had tried before your attempt!

I'm impressed with your knowledge of what operations are faster on the HP-41 and your use of synthetic register operations. It's been a while since I've done such programming, but it was lots of fun back when I had the time and energy to study it. When I get home, perhaps I can find my table of execution times for instructions.

The biggest speed improvement may be the backward testing, from high to low, as you mention above. We'll see if I can find time from work to tweak my code.


#58

Quote:
A very nice entry, and I like your ideas.
Thanks Steve for the kind words.
Quote:
Especially sweet as you apparently didn't look at what others had tried before your attempt!

Learned this from Valentine's challenges. There are thunderous intellects on this forum and so peeking ahead would rob the average Joe, I mean, Peter, from all the fun in solving things on my own.

Quote:
I'm impressed with your knowledge of what operations are faster on the HP-41 and your use of synthetic register operations. It's been a while since I've done such programming, but it was lots of fun back when I had the time and energy to study it. When I get home, perhaps I can find my table of execution times for instructions.

Don't be. The only calc I know something about is the 41 and I always liked speed improvements. You can actually still buy the Synthetix Pocket Reference guide which has execution time tables towards the end IIRC.

Quote:
The biggest speed improvement may be the backward testing, from high to low, as you mention above. We'll see if I can find time from work to tweak my code.

That would be great. I think if we combine the ideas that we have come up and incorporate it into the sleek algorithm from you/Hugh/Gerson we might have quite a speedy little program.

Cheers

Peter

Edited: 28 July 2009, 3:32 p.m.


#59

Timing for my RPL version above using an emulated 49G+ set to authentic speed is 18.7 seconds. I can chop that in half with one simple modification (using NEXTPRIME instead of INCR).

The reason why I was looking at improving the speed is to find the 7th number. So far, even after leaving the emu running at max speed for over 10 minutes, it hasn't been found yet. At max speed, the basic 6 number requirement takes 0.6 seconds so is about 29 times faster. That means that at real speed, finding the 7th number is going to take over 4h 50m !

Mark


#60

Is there a 7th number?

Not under 100,000,000 there isn't.


Pauli


#61

The Lucky Numbers of Euler
(by F. Le Lionnais)

Now that we have the programs, I started some digging to find out more about this sequence. I'm sure Egan will fill us in with far more clarity, but this is what I have learned so far.

Our little problem is related to Eulers Polynomial, even though Egan has disguised it a little. Note that the sequence x2-x gives the same result as x2+x just with the addition of 0 at the start, namely (0),2,6,12,20,30,...

Below is the quotation from wikipedia which talkes about a Heegener Number, something I have never heard off. Here is the link to Heegener numbers at Wolframs

It would appear that Rabinowitz proved that there is no result >41. The article is in german and I was not able to find a link to a free version on the internet.

Quote:
...from Wikipedia

Euler's prime-generating polynomial

n2 + n + 41,

which gives (distinct) primes for n=0,...,39, is related to the Heegner number 163 = 4*41 - 1.

Rabinowitz[1] proved that

n2 + n + p

gives primes for n=0,...,p-2 if and only if its discriminant 1-4p equals minus a Heegner number.

(Note that p-1 yields p2, so p-2 is maximal.) 1, 2, and 3 are not of the required form, so the Heegner numbers that work are 7,11,19,43,67,163, yielding prime generating functions of Euler's form for 2,3,5,11,17,41; these latter numbers are called lucky numbers of Euler by F. Le Lionnais.[2]





Edited: 29 July 2009, 11:04 a.m.


#62

It's to do with the class number of the imaginary quadratic field being 1. This is only the case for a few discriminants; Q[sqrt(-d)], for d in 3,4,7,8,11,19,43,67,163. if d = 1-4q, then we have q in; 2,3,5,11,17 and 41.

I originally wondered if there is a quick way to determine the class number h(1-4q) for small q. This would solve the problem, but it appears to need d = 1-4q to be factored (for the kronecker symbol see "class number" on mathworld).

if the class number is 2 we get some other polynomials

2x^2 + p, p in {3,5,11,29} all prime for x = 0,..p-1
2x^2+2x + n, n in {3,7,19} prime for x = 0,1,..n-2
px^2+px + n for pairs (p,n) in (3,2),(3,5),(3,11),(3,23),(5,3),(5,7),(5,13),(7,5),(7,17),(11,7),(13,11) prime for x = 0,1,..n-2

[from the "book of prime number records"]

#63

Nice response, thanks. I recognised the '41' as being Euler's polynomial but didn't remember the Euler bit and didn't notice the harmless change from x2 + x to x2 - x.

I guess my question is properly answered.

- Pauli

#64

Egan, this has been a very enjoyable challenge.

Just right for my level of skill and understanding, and of course the number 41 shows up prominently!

I appreciate the effort it takes to craft a challenge like this and present it to the community.

And it gave me an excuse to try out my recently acquired HP-41CX!

Thanks! And keep up the good work.

P.S. Oh, and good luck in the lottery too!

Edited: 29 July 2009, 5:56 p.m.


#65

Just want to echo Steve's words - thanks Egan, this has been most enjoyable. Especially due to the really fun code collaboration!

Thanks everyone, as usually, I learned a lot:-)

Cheers

Peter

#66

I know I was not allowed to do so but first I wrote a program in UserRPL:

\<< { } 1 ROT
FOR n
n DUP
IFERR
1 OVER 1 -
FOR i
IF DUP PRIM? NOT
THEN DOERR
END
i 2 * +
NEXT
THEN DROP
ELSE DROP +
END
NEXT
\>>

Run it with a number bigger than 41 and you will get:

{ 2 3 5 11 17 41 }

Then I used a Perl-script "rpl2focal.pl" to translate it into the following listing.
Nah, just kidding. Of course I had to use my brain. Without modules there's no such
thing as PRIM?. So we have to do that on our own. My goal was to make the innermost
loop as short as possible. That's why I keep the list of divisors in registers.
It might be faster though to keep only primes in that list but then it's not so easy
to calculate the upper bound in lines 46-51.

So here's my contribution. Now I'm going to read all the other solutions.

Egan, thanks a lot for the fun I had thanks to your challenge.

Best regards

Thomas Klemm


Listing

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

A Initialization

14
XEQ A

This program will fill the registers R00 - R14 with the following values:

R00:    2
R01: 3
R02: 5
R03: 7
R04: 11
R05: 13
R06: 17
R07: 19
R08: 23
R09: 25
R10: 29
R11: 31
R12: 35
R13: 37
R14: 41

You may have noticed that this is not exactly the list of primes.
It takes 0'02" to complete on an emulated HP-41CX.

B Lottery Numbers

14
XEQ B

This program will display the following numbers:

 5.0000
11.0000
17.0000
41.0000

It takes 1'45" to complete on an emulated HP-41CX.
So what about 2 and 3 you may ask and I must admit that's a bug.
But it's quiet easy to fix that: just add X#Y? after line 57.
After that it took 1'50" to complete.

C Is it a prime?

This program isn't really needed but you can use it to check
whether a number is a prime or not. Unless it is a prime 0 will
be returned, otherwise a value > 0.

Assuming the size is big enough you might run:

200
XEQ A

Now you can check numbers up to 358,801 = 599^2.
However it's probably not the fastest method.
It may work for bigger numbers too though for some it may not.

76543
XEQ C
91.0000

Hence it is a prime.

181663
XEQ C
0.0000

Hence it is not a prime. You can find the smallest factor in LAST x:

LASTX
389.0000

#67

A wonderful entry!

I like the way you used counts of prime candidates. I've been trying to put together something similar, skipping multiples of 3 as you did, but my efforts so far have taken 2m 50s or so.

I increased the number tested to 17 to get up to 49 as we did, since this is up to the assumed possible lottery numbers. It ran in 2m 10s at default emulator speed.

In my tests, turning the V41 emulator speed down to minimum most closely matched my genuine HP-41CX. It takes about twice as long to run at that setting on my computer at home. I'm not sure if that's exactly the same or not.

A very nice first attempt! I'm impressed.

P.S. My home computer is a bit slower than the one I was using earlier. Hard to calibrate these things precisely, I guess.

Edited: 1 Aug 2009, 11:42 p.m.

#68

I finally produced an acceptable version that stores an array to skip odd numbers divisible by 3.

It runs on my real HP-41CX in 2m 29.3s, but your mileage may vary.

As usual, improvements are welcome. Of course the challenge is officially over and my wife would like me to return to earth, so...

 01 LBL "41BD3"
02 CLX
03 SETSW
04 RUNSW
05 2
06 VIEW X
07 STO 02
08 XEQ 06
09 7.022
10 STO 05
11 3
12 STO 00
13 STO 03
14 LBL 01
15 SIGN
16 LBL 02
17 STO 01
18 X^2
19 LASTX
20 -
21 RCL 00
22 +
23 ENTER
24 ENTER
25 SQRT
26 RCL 03
27 STO 04
28 ST+ 04
29 LBL 03
30 X>Y?
31 GTO 05
32 R^
33 X<>Y
34 MOD
35 X=0?
36 GTO 09
37 SIGN
38 ST+ 04
39 CLX
40 RCL IND 04
41 GTO 03
42 LBL 05
43 SIGN
44 RCL 01
45 +
46 RCL 00
47 X<>Y
48 X#Y?
49 GTO 02
50 VIEW X
51 LBL 09
52 RCL IND 05
53 STO 00
54 ISG 05
55 GTO 01
56 STOPSW
57 BEEP
58 RCLSW
59 STOP
60 LBL 06
61 4
62 STO 03
63 22
64 49
65 LBL 07
66 STO IND Y
67 RCL 02
68 -
69 X<0?
70 RTN
71 DSE Y
72 STO IND Y
73 RCL 03
74 -
75 DSE Y
76 GTO 07
77 END
#69

Quotations

Quote:
although 49 is a given max for the problem here

Quote:
02 3.04902
03 STO 00

I wasn't aware of this fact. But I realized that there are no solutions bigger than 41.
That's why I considered that as an upper limit. From now on I will do my timing based on 49. IIRC the Swiss lottery used 6 out of 42.
Now they are using 45.


Quote:
  • as n is odd, x2-x + n is also always odd. So the primality test does not check for 2

That saved nearly 100 tests I guess. So my list is starting now with 3, 5, 7, ...


Quote:
04 2
05 VIEW X

I consider that a little cheating. But I'm sure I can use that trick as well since 2 and 3 are missing now in my solution.


Quote:
  • as digits are the slowest things on the 41, often used numbers are stored in registers

Quote:
  STO a  ;(3)

I wasn't aware that I could use register a without harm. Or probably I just forgot about that.
But since I'm not using subroutines I'm on the safe side. Using register P is probably a bad idea, I guess.


Quote:
Another optimisation would be to only test primes instead of odd numbers.

Not exactly what I did but a little closer. I'm using only 17 numbers instead of 24 of which only 25, 35 and 49 are composite.


Quote:
Lastly, I think that starting the x2-x+n test with n-1 is advantageous.
It is the largest number to be tested for prime and hence has the highest chance to be not prime, at which point no further tests have to be done.

I'm not sure whether I understood that correctly but if you suggest to count x down from n-1 to 0
then there are certain cases where you start testing a big prime but you could find the smallest
composite much faster.

E.g. n = 37

ascending:

37 is a prime, therefore you have to check 3 and 5 to find out

39 is composite, divisible by 3

3 checks needed to reject 37

descending:

1297 is a prime, therefore you have to check 3, 5, 7, ..., 35 to find out

1227 is composite, divisible by 3

18 checks needed to reject 37

Here are the other cases:

   7:  37
15: 197
21: 401
25: 577
27: 677

However as you may see from the following table there's only one case (namely 29) where you'd have to check more than on prime
prior to find the first composite. Primes are marked with a sharp #.

   3:    3#    5#
5: 5# 7# 11# 17#
7: 7# 9 13# 19# 27 37#
9: 9 11# 15 21 29# 39 51 65
11: 11# 13# 17# 23# 31# 41# 53# 67# 83# 101#
13: 13# 15 19# 25 33 43# 55 69 85 103# 123 145
15: 15 17# 21 27 35 45 57 71# 87 105 125 147 171 197#
17: 17# 19# 23# 29# 37# 47# 59# 73# 89# 107# 127# 149# 173# 199# 227# 257#
19: 19# 21 25 31# 39 49 61# 75 91 109# 129 151# 175 201 229# 259 291 325
21: 21 23# 27 33 41# 51 63 77 93 111 131# 153 177 203 231 261 293# 327 363 401#
23: 23# 25 29# 35 43# 53# 65 79# 95 113# 133 155 179# 205 233# 263# 295 329 365 403 443# 485
25: 25 27 31# 37# 45 55 67# 81 97# 115 135 157# 181# 207 235 265 297 331# 367# 405 445 487# 531 577#
27: 27 29# 33 39 47# 57 69 83# 99 117 137# 159 183 209 237 267 299 333 369 407 447 489 533 579 627 677#
29: 29# 31# 35 41# 49 59# 71# 85 101# 119 139# 161 185 211# 239# 269# 301 335 371 409# 449# 491# 535 581 629 679 731 785
31: 31# 33 37# 43# 51 61# 73# 87 103# 121 141 163# 187 213 241# 271# 303 337# 373# 411 451 493 537 583 631# 681 733# 787# 843 901
33: 33 35 39 45 53# 63 75 89# 105 123 143 165 189 215 243 273 305 339 375 413 453 495 539 585 633 683# 735 789 845 903 963 1025
35: 35 37# 41# 47# 55 65 77 91 107# 125 145 167# 191# 217 245 275 307# 341 377 415 455 497 541# 587# 635 685 737 791 847 905 965 1027 1091# 1157
37: 37# 39 43# 49 57 67# 79# 93 109# 127# 147 169 193# 219 247 277# 309 343 379# 417 457# 499# 543 589 637 687 739# 793 849 907# 967# 1029 1093# 1159 1227 1297#
39: 39 41# 45 51 59# 69 81 95 111 129 149# 171 195 221 249 279 311# 345 381 419# 459 501 545 591 639 689 741 795 851 909 969 1031# 1095 1161 1229# 1299 1371 1445
41: 41# 43# 47# 53# 61# 71# 83# 97# 113# 131# 151# 173# 197# 223# 251# 281# 313# 347# 383# 421# 461# 503# 547# 593# 641# 691# 743# 797# 853# 911# 971# 1033# 1097# 1163# 1231# 1301# 1373# 1447# 1523# 1601#
43: 43# 45 49 55 63 73# 85 99 115 133 153 175 199# 225 253 283# 315 349# 385 423 463# 505 549 595 643# 693 745 799 855 913 973 1035 1099 1165 1233 1303# 1375 1449 1525 1603 1683 1765
45: 45 47# 51 57 65 75 87 101# 117 135 155 177 201 227# 255 285 317# 351 387 425 465 507 551 597 645 695 747 801 857# 915 975 1037 1101 1167 1235 1305 1377 1451# 1527 1605 1685 1767 1851 1937
47: 47# 49 53# 59# 67# 77 89# 103# 119 137# 157# 179# 203 229# 257# 287 319 353# 389# 427 467# 509# 553 599# 647# 697 749 803 859# 917 977# 1039# 1103# 1169 1237# 1307# 1379 1453# 1529 1607# 1687 1769 1853 1939 2027# 2117
49: 49 51 55 61# 69 79# 91 105 121 139# 159 181# 205 231 259 289 321 355 391 429 469 511 555 601# 649 699 751# 805 861 919# 979 1041 1105 1171# 1239 1309 1381# 1455 1531# 1609# 1689 1771 1855 1941 2029# 2119 2211 2305

So I agree with you that chances are lower that bigger numbers are primes.
However if you run into one it takes much longer to prove it is one.

Listing

 01 LBL A                 24 LBL B                47 /                    
02 E3 25 E3 48 INT
03 / 26 STO a 49 RCL a
04 STO M 27 / 50 /
05 2 28 STO M 51 STO O
06 ENTER 29 LBL 02 52 LBL 04
07 X^2 30 RCL IND M 53 RDN
08 3 31 RCL X 54 RCL IND O
09 STO IND M 32 2 55 MOD
10 ISG M 33 - 56 X=0?
11 LASTX 34 RCL a 57 GTO 05
12 + 35 / 58 ISG O
13 GTO 00 36 STO N 59 GTO 04
14 LBL 01 37 RDN 60 RDN
15 RCL Z 38 LBL 03 61 ISG N
16 X<>Y 39 RCL N 62 GTO 03
17 R^ 40 INT 63 VIEW IND M
18 + 41 ST+ X 64 LBL 05
19 LBL 00 42 + 65 ISG M
20 STO IND M 43 ENTER 66 GTO 02
21 ISG M 44 ENTER 67 END
22 GTO 01 45 SQRT
23 RTN 46 3

Timing

16 XEQ A:     0'02.06"
16 XEQ B: 1'27.56"

Conclusions

I realized that my calculation of the upper bound in lines 45-48 could be improved:

 45 SQRT
46 1
47 -
48 3
49 /
50 INT

This would save around 30 tests in the innermost loop.
But since these lines are executed about 90 times the saving of time is nullified.

What I might try is to reverse all the indexes and use DSE instead of ISG.
This would save the division by E3. But then I'd have to solve the problem that 0 is skipped.

Just wanted to thank you for all the contributions so far. They have been a great help.
Nevertheless I don't think that I can beat Peter's solution. It seems that caching the
numbers isn't suitable in this situation. Most of the numbers can be rejected by simply
trying 3, 5 or 7 as divisors. But for exactly these cases my list doesn't differ from
what he's using.

Best regards

Thomas

#70

Only now I realized that I compared my program to Peter's buggy version from message #36. Compared to the fixed version in message #38 it isn't that bad. So I've concatenated my formerly separate parts A and B, added the missing view of 2 and 3, replaced subtraction of 2 with DSE, used register Q to store the number 3 and added the same stopwatch-code as in message #67 to measure the time.

Listing


01 LBL A 25 X<>Y 49 RCL Q
02 0 26 R^ 50 /
03 SETSW 27 + 51 INT
04 X<>Y 28 LBL 00 52 RCL a
05 RUNSW 29 STO IND O 53 /
06 E3 30 ISG O 54 STO O
07 STO a 31 GTO 01 55 LBL 04
08 / 32 LBL 02 56 RDN
09 STO M 33 RCL IND M 57 RCL IND O
10 STO O 34 ENTER 58 MOD
11 2 35 DSE X 59 X=0?
12 VIEW X 36 DSE X 60 GTO 05
13 ENTER 37 RCL a 61 ISG O
14 X^2 38 / 62 GTO 04
15 3 39 STO N 63 RDN
16 VIEW X 40 RDN 64 ISG N
17 STO Q 41 LBL 03 65 GTO 03
18 STO IND O 42 RCL N 66 VIEW IND M
19 ISG O 43 INT 67 LBL 05
20 LASTX 44 ST+ X 68 ISG M
21 + 45 + 69 GTO 02
22 GTO 00 46 ENTER 70 STOPSW
23 LBL 01 47 ENTER 71 RCLSW
24 RCL Z 48 SQRT 72 END

Timing

I used to pack the programs and ran them at least twice using the 2nd reading.

Peter's (#36+#38) : 01'38"
Steve's (#67) : 01'33"
mine (#69) : 01'28"
Egan's (#70) : 02'09"

Measured on a HP-41CX emulator. Maybe somebody could do that on a real calculator?


Edited: 2 Aug 2009, 10:08 p.m.


#71

I've been trying to benchmark and understand your code better. I'm having trouble in two spots. The main one is where you seem to be calculating the polynomial to test for primality:

 42 RCL N
43 INT
44 ST+ X
45 +

We've been testing for prime candidates P=x^2-x+N for given N. If I'm not mistaken you're testing P=x*2+N. When I made the change to your code, however, it did't seem to work. Probably I've missed something that should be obvious, or my programming changes were bad:

 42 RCL N
43 INT
44 X^2
45 LASTx
46 -
47 +

I'm also partly responsible for another possible problem area. In my suggested fix for Peter's entry, my code tested all x from 1 to N-2. As stated originally we should be testing from 0 (or 1) to N-1. This seems to be reflected in the code:

 35 DSE X
36 DSE X

I hope you can clarify these issues, or help me understand your ideas better. Thanks!

Steve


#72

Quote:
If I'm not mistaken you're testing P=x*2+N.

I'm calculating the sequence recursively, thus saving the multiplication:

dP(x) = P(x+1) - P(x) = 2x = x + x

What might be confusing is the fact that register N keeps track of x and not N. P(x) is kept only on the stack.

Quote:
As stated originally we should be testing from 0 (or 1) to N-1.

My code tests all x from 1 to N-1, however the counter in register N runs from 0 to N-2.

 x    dP(x)       P(x+1)
------------------------
0 0 0 + N
1 2 2 + N
2 4 6 + N
3 6 12 + N
4 8 20 + N
(...)

Hope this clarifies the issues.

Kind regards

Thomas


Edited: 4 Aug 2009, 8:20 a.m.

#73

It looks to me like there's a problem handling the stack in the loop at LBL 04.

I think you're trying to maintain the value P=x^2-x+N that we're checking for primality on the stack, but you've lost N that you need when you decide P is prime and loop to LBL 03 to check the next value of x.

I'm still trying to understand and see if there's an easy fix.


Edited: 3 Aug 2009, 2:21 p.m.


#74

Quote:
I think you're trying to maintain the value P=x^2-x+N that we're checking for primality on the stack

That's true.

Quote:
but you've lost N that you need when you decide P is prime and loop to LBL 03 to check the next value of x.

That's the other advantage of using recursion: N is kept hidden in P(x). Therefore no need to keep that value separately on the stack. Otherwise I couldn't take the full advantage of the most convenient feature that T is copied and not only moved to Z when executing MOD in line 58. The counter of the innermost loop (55-62) is kept outside the stack in register O thus not interfering. I can't think of a tighter solution.

Edited: 4 Aug 2009, 8:15 a.m.


#75

Very clever, Thomas!

Thanks for clearing that up.

I couldn't understand how you were getting the right answers without seeming to use the specified equation, but you've explained it well.

Sorry for being slow and assuming mistakes where there weren't any.

My "corrected version" isn't needed now, obviously.


#76

I had better commented my code as Gerson and Peter did. Most probably even I will have difficulties to understand my program in half a year.

#77

So what were my lottery numbers?

Steve was correct in assuming I selected this problem because it had a 41 in it. I also selected the lottery theme, as Peter and others have pointed out, because they are the only 6 known Euler Lucky Numbers (Fascinating!) and they happen to fall within most (if not all) lottery ranges. However they are not Egan's Lucky Numbers. That's right. I lost. Only one of the 6 (17) was lucky. I am a dollar poorer, but much wiser after reading this fascinating thread. Thank you all.

My solution was neither short, speedy, or sexy (read synthetic). I made no attempts to optimize after my initial proof-of-concept. I tip my party hat to all of you that participated with your brief, quick, and exotic solutions.

My solution was built on the following assumptions:

  1. If x2 - x + n is a prime number for x = 0, 1, ..., n - 1, then n must also be prime for the simple fact that when x = 0 or 1, then x2 - x = 0. So only prime n need be checked.
  2. Because of #1, there is only the need to check the range x = 2, 3, ..., n - 1.
  3. Because of #2, n = 2 can be a given (trivial), and it was given in my 2nd post. So you only need to find 5 odd prime numbers.

First I needed to generate a list of primes to be used as the candidates and for trial division. When I started on my prime number generator I started with the assumption that 2 was prime, I also knew that 3, 5, and 7 were also prime, but then so was 11. Before I knew it I had coded all 16 primes < 56 into my 41CX. It's amazing that being forced to memorize times-tables as a child that you were also indirectly memorizing a list of trivial primes. So I stuck with the bloated prime table method.

Now why < 56? I didn't give a range for the lottery because most range from 1-50 or greater. And the answer for any lottery known to me would have been the same. But, given that some searched from top down I should have given a range. I assumed that all would search from bottom up. So, back to why < 56? Well that was the range for the lottery I entered.

Once I had my prime table, my program simply processed each odd prime (#3), with the range in (#2). I used MOD instead of divide to check primality.

My program takes an embarrassingly 3.5 minutes on a 41CX with a printer attached (should run a bit faster without it). But the fun for me was watching others play.

BTW, I did dig a bit deeper into the mystery of Euler's Lucky Numbers, but did not find any quick or elegant tricks to mathematically make this problem easier to code. Such is life when working with prime numbers.

Again, thank you all for your participation.

Code:

 01 LBL "BDC2"       28 STO 10           55 RCL 18           82 ISG 18          
02 0 29 37 56 INT 83 GTO 01
03 SETSW 30 STO 11 57 X^2 84 RCL IND 17
04 FIX 00 31 41 58 LASTX 85 PRX
05 CF 29 32 STO 12 59 - 86 RCL 16
06 RUNSW 33 43 60 RCL IND 17 87 5
07 2 34 STO 13 61 + 88 X=Y?
08 STO 00 35 47 62 STO 19 89 GTO 05
09 3 36 STO 14 63 SQRT 90 1
10 STO 01 37 53 64 INT 91 ST+ 16
11 5 38 STO 15 65 STO 20 92 LBL 04
12 STO 02 39 2 66 0 93 ISG 17
13 7 40 PRX 67 STO 21 94 GTO 00
14 STO 03 41 1 68 LBL 02 95 LBL 05
15 11 42 STO 16 69 RCL 20 96 ADV
16 STO 04 43 1.01500 70 RCL IND 21 97 STOPSW
17 13 44 STO 17 71 X>Y? 98 FIX 06
18 STO 05 45 LBL 00 72 GTO 03 99 RCLSW
19 17 46 2 73 RCL 19 100 "TIME: "
20 STO 06 47 RCL IND 17 74 RCL IND 21 101 ATIME24
21 19 48 1 75 MOD 102 PRA
22 STO 07 49 - 76 X=0? 103 FIX 09
23 23 50 1000 77 GTO 04 104 SF 29
24 STO 08 51 / 78 1 105 END
25 29 52 + 79 ST+ 21
26 STO 09 53 STO 18 80 GTO 02
27 31 54 LBL 01 81 LBL 03

Edited: 2 Aug 2009, 7:50 p.m.


#78

Hi, Egan --

The losing Illinois Lottery Mega Millions ticket whose image you displayed -- and I doubt that you would have posted a winning one for the world to see -- seems to be one panel of the "pick six" Lotto Game type. 41 was your chosen "money ball" number.

According to the Illinois Lottery's website, the probability (not "odds", as is commonly misrepresented) of matching six unique numbers from 1 through 52 inclusive is 1:20,358,520.

Quote:

Lotto Game Panel - That area of an Official On-Line Lotto Game ticket identified by an alpha character, A through J, containing six (6) two-digit player or computer-selected numbers.

Lotto Game Winning Numbers - Six (6) two-digit numbers, from one (01) through fifty-two (52), randomly selected at each Lotto Game drawing, which shall be used to determine the winning Lotto plays contained on Official On-Line Lotto Game tickets.


http://www.illinoislottery.com/Games/OnLineRules.htm#_Toc123708150

20,358,520 is the statistical combination of 52 and 6.

The irony? Unfortunately, the HP-41 does not have statistical combination or permutation built-in, and the factorial function must be spelled out as "FACT". The HP-41 Stat Pac doesn't have 'em, either.


#79

Hello Karl,

IIRC, the Mega Millions is a multistate lottery with nearly 200 million combinations. The first 5 numbers are in the range 1-56 and the Mega Ball is in the range of 1-46. So, C(56,5)*46.

Quote:
The irony? Unfortunately, the HP-41 does not have statistical combination or permutation built-in, and the factorial function must be spelled out as "FACT". The HP-41 Stat Pac doesn't have 'em, either.

I commented on this in my last challenge. I was unable to find any HP41 ware that had COMB and building one from FACT was inaccurate for integer work. Quite a mystery to me that HP's flagship product for many years never had this function.

The good news is that HP41 development is not dead. SandMath 5 recently released with nCr and nPr.


#80

Angel Martin heard you with his new sandmath which has a mcode comb & perm!

cheers

peter

#81

This seems strange to me.

Looking at the lottery ticket printed, I get the impression that the payout is $77 million. But the probability of getting all 6 numbers correct is only 1 in 21 million!

That would seem like a good bet to make. However, it is well known that lottery payouts, on average, are terrible. Much worse than casino gambling for example. Usually the expected return is something like 50%. That's why it's such a money maker for whoever is funding the lottery (and those who take there cut along the way).

I'm sure I'm missing something. Can someone clear this up for me before I become a gambler?


#82

From above, its 1 in C(56,5)*46 or 194,810,616.


#83

Hmmh, my 15C as well as my 42S concur, displaying 175.711.536 = 1,76E8 instead.


#84

Oops. You are correct. Unsure where I got it wrong on Free42/iPhone. There is a risk to not having real keys.

Edit: Found it. It was an old calculation on the stack and I assumed was the value I precomputed earlier. C(56,6)*6 was the value assuming all six balls were 1-56 and that the power ball was the last drawn. 1/C(56,6) is the probability of getting all 6 balls correct with a 1/6 chance of selecting the correct ball as the power ball. C(56,5)*51 is the same.

Edited: 3 Aug 2009, 1:53 p.m.

#85

194,810,616 is C(56,5)*51

Mark


#86

The power ball is a duplicate ball in the range 1-46. Something I failed to report earlier in my summary. I'm clearly a novice at this lottery stuff. I live in a state with no lottery (Utah). When I staged the challenge I knew I was going to be in Chicago that week. So, I made a mental note to get a ticket for a lottery ending on the 31st as I said I would in the challenge.

There were many lotteries to choose from, but I liked the size of the pot for the Mega Millions and it also allowed me to choose 41 as a position of prominence.

So the combinations are C(56,5) for the main 5 balls and 46 for the mega ball. With a total of C(56,5)*46. Apparently the two sets can intersect as seen if you Google for Mega Millions past numbers.

Edited: 3 Aug 2009, 2:06 p.m.

#87

At $77M Mega Millions is a mega-bad bet, but there are times in most lottery games where the bet becomes favorable in that your expected payout is greater is greater than your bet. The problem is that lotteries don't publish information on their current take so you can't accurately figure this out. You know the prize and you know the odds of winning but you don't know your odds of having to share that prize.

Still, I often join the throngs of idiots playing the lottery in the hope of winning enough money to acquire an HP-95C some day :)


#88

I buy a ticket now and then as well. I figure that buying a ticket increases my odds of winning, but not by much. There's a still a much greater chance of getting struck by lightning, which is still free. I'd rather win the lottery, so I figure a ticket improves my odds relative to either winning the lottery or being struck by lightning....

#89

All lotteries are a bad bet. Even when the pot gets big so that it starts to seem like a good idea, the reason it is getting that big is there are a lot of tickets purchased! Therefore your odds of having to share the pot increase proportionally.

Those who run the system are not stupid. They don't let the odds get in the players favor, even when they make it seem like it is.

I just wanted to make sure everyone is clear on that. It can be fun to play, but you shouldn't ever bet a lot assuming your expectation of winning is good enough to justify the expense. Fund your IRA or 401K or make other sound investments instead. It's what anyone who understands math does.


#90

Quote:
Those who run the system are not stupid. They don't let the odds get in the players favor, even when they make it seem like it is.

I agree with your advice that the money should simply go to savings/investments. I've never played a lottery.

However, from an analysis point of view the expected value for your $1 lotto ticket does increase as the jackpot rolls over, and the interesting question is when is the expected value at least $1 for your purchased ticket.

The lotto people are not stupid, as you say, they get a fixed amount of money from the sale of tickets, usually just under 50%, this is true even with the fixed payout prizes since the fixed payouts are deducted from the jackpot. So regardless of increasing player expected value, the lotto people get their money always, and consistently!

Exploring the point at which the expected value reaches at least $1 would take into account adjusting the jackpot amount with:

 * Roll over jackpot amount
* Number of tickets purchased on current go around.
* Chance the jackpot needs to be split (derived from above)
* During the current game how much money goes toward the jackpot with each ticket (Probably about 48%).
* The Jackpot cash payout.. not the annuity (the cash payout is about half the reported jackpot)
* Taxes

My feeling is that the point at which the expected value is equal or greater to the purchase price for a ticket is much higher then people think. Something for me to work on on a rainy day :)


#91

How about "fantasized value" or "perceived payoff"?

Being told of an even larger pot of gold at the end of the rainbow doesn't make the search any more likely to succeed. (It's an imperfect analogy, but not by much!)


#92

Quote:
How about "fantasized value" or "perceived payoff"?

http://en.wikipedia.org/wiki/Expected_value

Edited: 9 Aug 2009, 7:09 p.m.

#93

Quote:
Still, I often join the throngs of idiots playing the lottery in the hope of winning enough money to acquire an HP-95C some day :)

Isn't that what an eBay watch list is for? :-)
#94

I view lotteries as a tax on the innumerate and lots of people participating is good -- I pay less other tax.

- Pauli

#95

Inspired by your table of primes I wondered whether I could use that to make my program faster.

The biggest problem was to determine the upper bound for the innermost loop (73-78). This might be improved since there are still cycles that could be avoided.

Using + to create the primes turned out to be nearly a second faster than just entering them.

Listing

 01 LBL "BD"                   26 +                          51 +                          76 RCL 00   ; 3         
02 0 27 STO 03 ; 11 52 STO 12 ; 43 77 /
03 SETSW 28 + 53 + 78 INT
04 X<>Y 29 STO 04 ; 13 54 + 79 RCL 15 ; 1000
05 RUNSW 30 + 55 STO 13 ; 47 80 /
06 E3 31 + 56 LBL 02 81 STO O
07 STO 15 ; 1000 32 STO 05 ; 17 57 RCL IND M 82 LBL 04
08 / 33 + 58 RCL a ; 2 83 RDN
09 STO M 34 STO 06 ; 19 59 - 84 RCL IND O
10 SIGN 35 + 60 RCL 15 ; 1000 85 INT
11 STO 14 ; 1 36 + 61 / 86 MOD
12 ENTER 37 STO 07 ; 23 62 RCL 14 ; 1 87 X=0?
13 + 38 + 63 + 88 GTO 05
14 ENTER 39 + 64 STO N 89 ISG O
15 ENTER 40 + 65 RCL IND M 90 GTO 04
16 STO a ; 2 41 STO 08 ; 29 66 LBL 03 91 RDN
17 VIEW X 42 + 67 RCL N 92 ISG N
18 LASTX 43 STO 09 ; 31 68 INT 93 GTO 03
19 + 44 + 69 ST+ X 94 VIEW IND M
20 STO 00 ; 3 45 + 70 + 95 LBL 05
21 + 46 + 71 ENTER 96 ISG M
22 STO 01 ; 5 47 STO 10 ; 37 72 ENTER 97 GTO 02
23 + 48 + 73 SQRT 98 STOPSW
24 STO 02 ; 7 49 + 74 RCL 14 ; 1 99 RCLSW
25 + 50 STO 11 ; 41 75 - 100 END

13
XEQ "BD"

Timimg: 01'22.19"


#96

Quote:
Using + to create the primes turned out to be nearly a second faster than just entering them.

Interesting that there is a penalty when storing constants in 41 RPN code.

#97

As Peter stated:

Quote:
Digits really are very slow on the hp41.

#98

Just a curiosity: lines 07 through 22 can be replaced with the code below, thus saving three steps.

07 .00701
08 STO 17
09 2.566543832
10 LBL 05
11 ENTER
12 INT
13 STO IND 17
14 /
15 1
16 -
17 1/X
18 ISG 17
19 GTO 05


Edited: 8 Aug 2009, 12:19 p.m.


#99

After what Peter and Thomas have stated about numbers in 41 RPN code, I expect that it will run faster too.

I have another challenge brewing and it will require many more prime numbers. I guess I'll have to write one after all.


Quote:
After what Peter and Thomas have stated about numbers in 41 RPN code, I expect that it will run faster too.

I don't believe so because of the constant in line 09. I remember Peter told something about number length too. This makes sense because on the 15C, for instance, that constant alone would take up eleven steps. Anyway, this is just an exotic way to generate a (somewhat limited) prime table (http://hypography.com/forums/physics-and-mathematics/17006-the-holy-grail-of-mathematics-8.html).

In fact, your program could be seven steps shorter, but it would require a little more time to run, as we can see by the output:


01 LBL "BD" 28 2 37 RCL 18 64 ISG 18
02 0 29 RCL IND 17 38 INT 65 GTO 01
03 SETSW 30 1 39 X^2 66 RCL IND 17
04 FIX 0 31 - 40 LASTX 67 PRX
05 CF 29 32 1000 41 - 68 RCL 16
06 RUNSW 33 / 42 RCL IND 17 69 5
07 .00701 34 + 43 + 70 X=Y?
08 STO 17 35 STO 18 44 STO 19 71 GTO 05
09 2.566543832 36 LBL 01 45 SQRT 72 1
10 XEQ "PT" ------------- 46 INT 73 ST+ 16
11 1.36280923 47 STO 20 74 LBL 04
12 + 48 0 75 ISG 17
13 .006 49 STO 21 76 GTO 00
14 ST+ 17 50 LBL 02 77 LBL 05
15 x<>y 51 RCL 20 78 ADV
16 XEQ "PT" ------------- 52 RCL IND 21 79 STOPSW
17 47 01 LBL "PT" 53 X>Y? 80 FIX 6
18 STO 14 02 ENTER 54 GTO 03 81 RCLSW
19 53 03 INT 55 RCL 19 82 "TIME: "
20 STO 15 04 STO IND 17 56 RCL IND 21 83 ATIME24
21 2 05 / 57 MOD 84 PRA
22 PRX 06 1 58 X=0? 85 FIX 9
23 1 07 - 59 GTO 04 86 SF 29
24 STO 16 08 1/X 60 1 87 END
25 1.01500 09 ISG 17 61 ST+ 21
26 STO 17 10 GTO "PT" 62 GTO 02
27 LBL 00 11 END 63 LBL 03

------------------------

2 ***
3 ***
5 ***
11 ***
17 ***
41 ***

TIME: 00:04:08.89

------------------------


Edited: 8 Aug 2009, 4:41 p.m.


Possibly Related Threads...
Thread Author Replies Views Last Post
  Fun graphs on HP Prime Mic 9 1,110 09-15-2013, 08:30 AM
Last Post: Eddie W. Shore
  Fun things found by running strings on the 39gII emulator bhtooefr 11 1,568 05-16-2013, 12:40 AM
Last Post: Mic
  My birthday, so a little commemorative mini-challenge ! Valentin Albillo 43 3,907 03-07-2013, 03:44 AM
Last Post: Walter B
  HP-28S plays "Happy Birthday" Kevin F 0 398 08-06-2012, 04:53 PM
Last Post: Kevin F
  Fun with the HP-35 Gerson W. Barbosa 11 1,361 08-02-2012, 10:58 AM
Last Post: Gerson W. Barbosa
  First post here: something fun... Jeff Dinkins 4 683 04-22-2012, 02:42 PM
Last Post: Jeroen Van Nieuwenhove
  devolution of the wp34s label---or fun with colors troy 19 2,311 02-09-2012, 02:52 PM
Last Post: Egan Ford
  Happy 40th Birthday, HP35 Jake Schwartz 8 1,082 01-05-2012, 10:42 PM
Last Post: db (martinez, ca.)
  HP 15c searchable Adv Fun Handbook (link) - done by Joe Horn Gene Wright 19 1,855 11-24-2011, 04:59 AM
Last Post: Joel Setton (France)
  Happy Birthday HP9810A Achim (Germany) 6 855 11-03-2011, 03:12 PM
Last Post: Achim (Germany)

Forum Jump: