MathGroup Archive 2007

[Date Index] [Thread Index] [Author Index]

Search the Archive

Re: compile speed

  • To: mathgroup at smc.vnet.net
  • Subject: [mg74530] Re: compile speed
  • From: "Boson" <sandro.romani at gmail.com>
  • Date: Sat, 24 Mar 2007 05:21:15 -0500 (EST)
  • References: <etqnm3$nf$1@smc.vnet.net><ett7ba$hll$1@smc.vnet.net>

On Mar 22, 7:22 am, "Ray Koopman" <koop... at sfu.ca> wrote:
> On Mar 21, 12:42 am, "Boson" <sandro.rom... at gmail.com> wrote:
>
>
>
> > dear mathematica users,
>
> > i've written a simple function that works on a pair of binary
> > matrices:
> > (mathematica 5.2 linux, on a 32 bit platform)
>
> > tab[nx_, ny_, frac_] := Table[If[Random[] < frac, 1, 0], {nx}, {ny}]
>
> > nx = 25; ny = 50; frac1 = 0.1; frac2 = 0.5;
> > p1 = 0.4; p2 = 0.2;
> > tabrect = tab[nx, ny, frac1];
> > tabsq = tab[ny, ny, frac2];
>
> > testnocomp[mat1_, mat2_, n1_, n2_, pp1_, pp2_] := Module[{tmp, sum,
> > val},
> >     tmp = mat2; Do[sum = mat1[[k,j]] + mat2[[k,i]];
> >        val = Which[sum == 2, If[Random[] < pp1, 1, tmp[[i,j]]], sum ==
> > 1,
> >          If[Random[] < pp2, 0, tmp[[i,j]]], sum == 0, tmp[[i,j]]];
> > tmp[[i,j]] = val,
> >       {k, n1}, {i, n2}, {j, n2}]; tmp];
>
> > Timing[resnc = testnocomp[tabrect, tabsq, nx, ny, p1, p2]; ]
>
> > the result of the timing is
> > {0.7558840000000013*Second, Null}
>
> > since i need high values of nx,ny (~5000) and the loop scales as
> > nx*ny^2,
> > i tried to implement a compiled version of the previous function:
>
> > test := Compile[{{mat1, _Integer, 2}, {mat2, _Integer, 2}, {n1,
> > _Integer}, {n2, _Integer},
> >      {pp1, _Real}, {pp2, _Real}}, Module[{tmp, sum, val},
> >      tmp = mat2; Do[sum = mat1[[k,j]] + mat2[[k,i]];
> >         val = Which[sum == 2, If[Random[] < pp1, 1, tmp[[i,j]]], sum
> > == 1,
> >           If[Random[] < pp2, 0, tmp[[i,j]]], sum == 0, tmp[[i,j]]];
> > tmp[[i,j]] = val,
> >        {k, n1}, {i, n2}, {j, n2}]; tmp], {{Random[_], _Real}}];
>
> > Timing[res = test[tabrect, tabsq, nx, ny, p1, p2]; ]
>
> > results are a disaster:
> > {14.814747999999996*Second, Null}
>
> > i'm sure this is related to my poor mathematica programming
> > experience..
>
> > could you suggest me a faster version to solve this problem?
>
> > regards,
> > sandro
>
> First, whenever you Compile a function that is at all complicated
> (and perhaps always), you should look at its InputForm:
>
> InputForm[test]
>
> CompiledFunction[{{_Integer, 2}, {_Integer, 2}, _Integer,
>   _Integer, _Real, _Real}, {{2, 2, 0}, {2, 2, 1}, {2, 0, 0},
>   {2, 0, 1}, {3, 0, 0}, {3, 0, 1}, {2, 2, 2}},
>  {0, 10, 3, 0, 3}, {{1, 5}, {12, 1, 2}, {9, 0, 2}, {9, 1, 3},
>   {9, 1, 4}, {4, 0, 5}, {79, 5, 2, 11}, {4, 0, 6},
>   {79, 6, 3, -2}, {4, 0, 7}, {79, 7, 4, -2},
>   {64, 0, 0, 5, 0, 7, 0, 8}, {64, 1, 0, 5, 0, 6, 0, 9},
>   {24, 8, 9, 8}, {21, Function[{mat1, mat2, n1, n2, pp1,
>      pp2}, Which[sum == 2, If[Random[] < pp1, 1, tmp[[i,j]]],
>      sum == 1, If[Random[] < pp2, 0, tmp[[i,j]]], sum == 0,
>      tmp[[i,j]]]], {sum, 2, 0, 8, Module},
>    {i, 2, 0, 6, Block}, {tmp, 2, 2, 2, Module},
>    {j, 2, 0, 7, Block}, {k, 2, 0, 5, Block}, 2, 2, 0, 2, 2,
>    1, 2, 0, 0, 2, 0, 1, 3, 0, 0, 3, 0, 1, 3, 0, 2},
>   {21, Function[{mat1, mat2, n1, n2, pp1, pp2},
>     tmp[[i,j]] = val], {val, 3, 0, 2, Module},
>    {i, 2, 0, 6, Block}, {tmp, 2, 2, 2, Module},
>    {j, 2, 0, 7, Block}, {k, 2, 0, 5, Block}, 2, 2, 0, 2, 2,
>    1, 2, 0, 0, 2, 0, 1, 3, 0, 0, 3, 0, 1, 6, 0, 17},
>   {42, -6}, {2}}, Function[{mat1, mat2, n1, n2, pp1, pp2},
>   Module[{tmp, sum, val}, tmp = mat2;
>     Do[sum = mat1[[k,j]] + mat2[[k,i]];
>       val = Which[sum == 2, If[Random[] < pp1, 1,
>          tmp[[i,j]]], sum == 1, If[Random[] < pp2, 0,
>          tmp[[i,j]]], sum == 0, tmp[[i,j]]];
>       tmp[[i,j]] = val, {k, n1}, {i, n2}, {j, n2}]; tmp]],
>  Evaluate]
>
> If there is ordinary Mathematica code in the middle of the compiled
> code (other than the block at the end, which is always there), then
> the compiled code will probably be slower than the uncompiled code.
>
> Here are three successively faster versions of the uncompiled code.
> All three put the k-loop inside the i- and j-loops and give the same
> result, which necessarily differs from yours because changing the
> loops changes the random numbers. However, the new versions are
> statistically equivalent to yours.
>
> tab[nx_, ny_, frac_] := Table[If[Random[] < frac, 1, 0], {nx}, {ny}]
>
> nx = 25; ny = 50; frac1 = 0.1; frac2 = 0.5;
> p1 = 0.4; p2 = 0.2;
> tabrect = tab[nx, ny, frac1];
> tabsq = tab[ny, ny, frac2];
>
> testnocomp[mat1_, mat2_, n1_, n2_, pp1_, pp2_] := Module[
> {tmp, sum, val},
>     tmp = mat2; Do[sum = mat1[[k,j]] + mat2[[k,i]];
>        val = Which[sum == 2, If[Random[] < pp1, 1, tmp[[i,j]]],
>                    sum == 1, If[Random[] < pp2, 0, tmp[[i,j]]],
>                    sum == 0, tmp[[i,j]]];
>        tmp[[i,j]] = val,
>       {k, n1}, {i, n2}, {j, n2}]; tmp];
>
> Timing[resnc = testnocomp[tabrect, tabsq, nx, ny, p1, p2]; ]
> {1.26 Second,Null}
>
> testnocomp2[mat1_, mat2_, n1_, n2_, pp1_, pp2_] := Module[{tmp, sum},
> Table[sum = mat1[[All,j]] + mat2[[Range@n1,i]]; tmp = mat2[[i,j]];
>       Do[Which[sum[[k]] == 2, If[Random[] < pp1, tmp = 1],
>                sum[[k]] == 1, If[Random[] < pp2, tmp = 0]],{k, n1}];
>       tmp, {i, n2}, {j, n2}]];
>
> SeedRandom[1];
> Timing[resnc2 = testnocomp2[tabrect, tabsq, nx, ny, p1, p2]; ]
> {0.67 Second,Null}
>
> testnocomp3[mat1_, mat2_, n1_, n2_, pp1_, pp2_] := Module[{tmp},
> Table[tmp = mat2[[i,j]];
>       Scan[Which[# == 2, If[Random[] < pp1, tmp = 1],
>                  # == 1, If[Random[] < pp2, tmp = 0]]&,
>            mat1[[All,j]] + mat2[[Range@n1,i]]];
>       tmp, {i, n2}, {j, n2}]];
>
> SeedRandom[1];
> Timing[resnc3 = testnocomp3[tabrect, tabsq, nx, ny, p1, p2]; ]
> {0.57 Second,Null}
>
> testnocomp4[mat1_, mat2_, n1_, n2_, pp1_, pp2_] :=
> Table[Fold[Which[#2 == 2, If[Random[] < pp1, 1, #1],
>                  #2 == 1, If[Random[] < pp2, 0, #1],
>                  True, #1]&,
>            mat2[[i,j]], mat1[[All,j]] + mat2[[Range@n1,i]]],
>       {i, n2}, {j, n2}];
>
> SeedRandom[1];
> Timing[resnc4 = testnocomp4[tabrect, tabsq, nx, ny, p1, p2]; ]
> {0.09 Second,Null}
>
> resnc2 === resnc3 === resnc4
> True

thanks everybody!
the fastest among suggestions is only a factor of 10 slower than a
compiled C version.
i'd say this is not a bad result at all, but i'll keep experimenting.

i'm studying the code you all posted
(functional programming and #& are really powerful, but i still need
to absorb them)

thanks again,
Sandro



  • Prev by Date: Re: Which Mathematica product should I get?
  • Next by Date: Some questions about vector
  • Previous by thread: Re: compile speed
  • Next by thread: rescale gray levels with ListDensityPlot