Question: Playing Poker in Maple

I have this code (see below) which calculates the probability of finding 4 aces in a 52 card deck through simulation.

My problem is however that the code is very slow and the probability is very low which means that I have to run the

simulation maybe one milion times to get an good approximation. The code can not handle this.

Is there any way to make the code faster so it can handle one million simulations ?

 

restart:
randomize():
with(combinat, randcomb):
with(ListTools):

Clubs := [C[A], C[K], C[Q], C[J], seq(C[i], i = 2 .. 10)] :
Diamonds := [D[A], D[K], D[Q], D[J], seq(D[i], i = 2 .. 10)] :
Hearts := [H[A], H[K], H[Q], H[J], seq(H[i], i = 2 .. 10)] :
Spades := [S[A], S[K], S[Q], S[J], seq(S[i], i = 2 .. 10)] :
Cards := [op(Clubs), op(Diamonds), op(Hearts), op(Spades)] :

n := 1000 :

x1 := [seq(randcomb(Cards, 4), i = 1 .. n)] :

for i to n do

if 

Occurrences(C[A], x1[i]) = 1  and 
Occurrences(D[A], x1[i]) = 1  and
Occurrences(H[A], x1[i]) = 1  and 
Occurrences(S[A], x1[i]) = 1  then  CC[i] := 1:  else CC[i] := 0 :

end if end do :

"Probability 4 Aces" = evalf(add(CC[i], i = 1 .. n)/n, 8) ;

 

Please Wait...