# Code efficiency improvement needed solving problem 499 in Project Euler

I am trying to solve https://projecteuler.net/problem=499 by constructing a transition matrix and finding its steady state. Here is my code:

pmatrix[s_, m_, size_, iteration_] := Module[{matrix, matrix2},
matrixtransformer[n_] := 0;
Do[matrixtransformer[2^(i - 1) - m] = N[1/2^i,20], {i,
Floor[Log[2, size]]}];
matrixtransformer2[n_] := 0; matrixtransformer2[0] = 1;
matrix = Join[
Map[matrixtransformer2,
Table[Table[i - j, {i, size}], {j, m - 1}], {2}],
Map[matrixtransformer,
Table[Table[i - j, {i, size}], {j, m, size}], {2}]];
(*construct the initial matrix*)
fastmul[n_, t_] :=
Module[{matrix3}, matrix3 = n;
Do[matrix3 = Dot[matrix3, matrix3], {t}]; matrix3];
matrix = fastmul[matrix, iteration];
matrix2 =
Map[matrixtransformer2, Table[s - i, {i, 1, size}]].matrix;
N[1 - Plus @@ matrix2[[1 ;; m - 1]], 20]]


It works just fine with small parameters like $pmatrix[2,2,100,500]=0.25220849$ However, it is impossible to do the problem with matrix size of $10^9$. What can be done to improve the code? (PS:I am new with Wolfram Language, any suggestions on optimising the code will be appreciated also)

Update: I wrote a new algorithm seems to have a much better time complexity. But the ten thousand case would take a long time also, I don't know where the problem is. Here is my code:

pmatrix2[s_, m_, acc1_, acc2_, iteration_] :=
Module[{}, ptemp[n_] := 0;
ft[n_] := Sum[1/2^i ptemp[n + 2^(i - 1) - m], {i, acc2}];
Do[ptemp[i] = 1, {i, m - 1}];
Do[Do[ptemp[i] = ft[i], {i, m, acc1}], {iteration}];
1 - ptemp[s] // N]

• @Mr.Wizard it is up to 501 now. – happy fish Feb 6 '15 at 14:02
• Could someone please explain the sentence " The gambler is certain to win at least 1 pound, the starting value of the pot, at the cost of m pounds, the initial fee." to me? (I feel too dumb or too drunk, undecided yet) – Dr. belisarius Feb 6 '15 at 14:32
• @belisarius That means, when you spend m pounds to play the game, the worst case is that you have tail at the first time, then you collect the initial value of the pot, 1 pounds. That is "win at least 1 pound at the cost of m pounds". If you have head instead, the value will double until you have tail, then you collect the value of that time. – happy fish Feb 6 '15 at 14:52
• Oh, thanks. I find the writing quite obscure. – Dr. belisarius Feb 6 '15 at 15:12
• How is p2(2) roughly .2522? I pay my full stake on the first game. Heads I get nothing, tails I get 1 and game over. Either way I'm done. – Daniel Lichtblau Feb 6 '15 at 19:06

Making some suggestions for you code...

Do[matrixtransformer[2^(i - 1) - m] = N[1/2^i,20], {i,
Floor[Log[2, size]]}];


Note that I've used N[] for the values. Since we are not interested in the exact values using N[] reduces time.

And for the second part, you can use MatrixPower[] instead of iterating the multiplication. Its quiet fast. But the big problem is not there... That is with the size of the matrix.. Am also stuck here... But I think the trick is figuring out the Steady state matrix by using more math than programming... Anyways, the problem is enriching in so many ways.. so quiet happy to see a code here...

• Hi ! As it is currently written it is not a full-blown answer rather than a comment. If you, for example, make a side-by-side pre/after and post it with timings that would be great ! – Sektor Feb 6 '15 at 21:24
• MatrixPower[] reduces AbsoluteTiming by a factor of six for pmatrix[2, 2, 100, 500]. – bbgodfrey Feb 7 '15 at 2:57
• Thank you, I have taken your advice and updated the code – happy fish Feb 7 '15 at 4:59
• Run times for parameters 2, 2, 100, 500: - pmatrix as written in Question 1.65 sec - pmatrix with MatrixPower 0.27 sec - pmatrix with MatrixPower and N gave worse run times, even with precision less than \$MachinePrecision What are the values of acc1 and acc2 for your test problem? – bbgodfrey Feb 7 '15 at 5:23