Significant slowdown in a generalised implementation

Question

Consider the following code:

c0 := -1.24; 
b := 0.125; 
hh[{x_, y_}] := {x^2 + c0 - b*y, x}; 
hh2[{x_, y_}] := hh[hh[{x, y}]]; 
dhh := Transpose[{D[hh[{x, y}], x], D[hh[{x, y}], y]}]; 
dhh2 := Transpose[{D[hh[hh[{x, y}]], x], D[hh[hh[{x, y}]], y]}]; 
fixpoints = NSolve[hh[{x, y}] == {x, y}, {x, y}]; 
period2cycle = NSolve[hh[hh[{x, y}]] == {x, y}, {x, y}]; 
center = {x, y} /. fixpoints[[1]]; 
centerOther = {x, y} /. fixpoints[[2]]; 
pushvector = Eigensystem[dhh /. fixpoints[[1]]][[2,1]]; 
pushvectorOther = Eigensystem[dhh /. fixpoints[[2]]][[2,1]]; 
sink = {x, y} /. period2cycle[[1]]; 
pushvector2 = Eigensystem[dhh /. fixpoints[[2]]][[2,1]]; 
maxIter = 63; 
maxXXSize = 600; 
maxYYSize = 600; 
smalldelta = 0.05/maxXXSize; 
newfun[{{x_, y_}, n_}] := {hh[{x, y}], n + 1}; 
newtest[{{x_, y_}, n_}] := Abs[x] < 10000 && 
    Norm[{x, y} - sink] > 0.001; 
newcolorassign[{{x_, y_}, n_}] := Which[Abs[x] < 3, {0, 0, 0}, 
    Mod[Floor[Arg[x]/2^16], 2] > 0, {0.5, 0.5, 0.5}, True, 
    {1, 1, 1}];

With this setup one can execute the following:

Graphics[{Raster[Table[newcolorassign[NestWhile[newfun, 
       {center - pushvector*smalldelta*I*(k - maxXXSize/2 + 
           I*(j - maxYYSize/2)), 0}, newtest, 1, maxIter]], 
     {k, 1, maxXXSize}, {j, 1, maxYYSize}]], 
   PointSize[12/maxXXSize], White, 
   Point[{(maxXXSize - 1)/2, (maxYYSize - 1)/2}], 
   PointSize[4/maxXXSize], Black, 
   Point[{(maxXXSize - 1)/2, (maxYYSize - 1)/2}]}]

In a notebook this will produce a picture but because you specify the parameters at the beginning, to avoid confusion you should start a new notebook for each parameter pair (b,c0). Now to get around this problem I tried to generalise the code to the following.

hh1[{x_, y_}, b_, c_] := {x^2 + c - b*y, x}; 
hhhh[{x_, y_}, b_, c_] := hh1[hh1[{x, y}, b, c], b, c]; 
dhhh[b_, c_] := Transpose[{D[hh1[{x, y}, b, c], x], 
     D[hh1[{x, y}, b, c], y]}]; 
dhhh2[b_, c_] := Transpose[{D[hhhh[{x, y}, b, c], x], 
     D[hhhh[{x, y}, b, c], y]}]; 
fixpts[b_, c_] := NSolve[hh1[{x, y}, b, c] == {x, y}, {x, y}]
per2cycle[b_, c_] := NSolve[hhhh[{x, y}, b, c] == {x, y}, {x, y}]
centre[b_, c_] := {x, y} /. fixpts[b, c][[1]]; 
centreOther[b_, c_] := {x, y} /. fixpts[b, c][[2]]; 
pushvect[b_, c_] := Eigensystem[dhhh[b, c] /. 
      fixpts[b, c][[1]]][[2,1]]; 
pushvectOther[b_, c_] := 
   Eigensystem[dhhh[b, c] /. fixpts[b, c][[2]]][[2,1]]; 
snk[b_, c_] := {x, y} /. per2cycle[b, c][[1]]; 
pushvect2[b_, c_] := Eigensystem[dhhh[b, c] /. 
      fixpts[b, c][[2]]][[2,1]]; 
maxIter = 63; 
maxXXSize = 600; 
maxYYSize = 600; 
smalldelta = 0.05/maxXXSize; 
newfun[{{x_, y_}, n_}, b_, c_] := {hh1[{x, y}, b, c], n + 1}; 
newtest[{{x_, y_}, n_}, b_, c_] := Abs[x] < 10000 && 
    Norm[{x, y} - snk[b, c]] > 0.001; 
newcolorassign[{{x_, y_}, n_}] := Which[Abs[x] < 3, {0, 0, 0}, 
    Mod[Floor[Arg[x]/2^16], 2] > 0, {0.5, 0.5, 0.5}, True, 
    {1, 1, 1}]; 
unstablepic[b_, c_] := Graphics[
    {Raster[Table[newcolorassign[NestWhile[newfun[#1, b, c] & , 
         {centre[b, c] - pushvect[b, c]*smalldelta*I*
            (k - maxXXSize/2 + I*(j - maxYYSize/2)), 0}, 
         newtest[#1, b, c] & , 1, maxIter]], {k, 1, maxXXSize}, 
       {j, 1, maxYYSize}]], PointSize[12/maxXXSize], White, 
     Point[{(maxXXSize - 1)/2, (maxYYSize - 1)/2}], 
     PointSize[4/maxXXSize], Black, 
     Point[{(maxXXSize - 1)/2, (maxYYSize - 1)/2}]}];

Now unstablepic[0.125,-1.24] produces the same output as the original code but is significantly slower. To see this more quickly, one can take maxXXSize = 60 and maxYYSize = 60 instead of 600. I’m not sure why the generalised code is taking so much longer than the original. Is there a way to make this quicker?

performance tuning

b3m2a1 · Accepted Answer

It looks to me like you're recalculating a lot of stuff that should be constant (unless you're doing dangerous things with global variables).
Here's a way to not do that, mostly by changing newTest and
newtest[{{x_, y_}, n_}, b_, c_, snk_] :=
  Abs[x] < 10000 && Norm[{x, y} - snk] > 0.001;
unstablePts[b_, c_] :=
 With[
  {
   centre = centre[b, c], 
   pushvect = pushvect[b, c],
   sink = snk[b, c]
   },
  Table[
   NestWhile[
    newfun[#1, b, c] &,
    {
     centre - pushvect*
       smalldelta*I*(k - maxXXSize/2 + I*(j - maxYYSize/2)), 0},
    newtest[#1, b, c, sink] &,
    1,
    maxIter
    ],
   {k, 1, maxXXSize},
   {j, 1, maxYYSize}
   ]
  ]

unstablepic[b_, c_] :=
  Graphics[{
    Raster[Map[newcolorassign, unstablePts[b, c], {2}]],
    PointSize[12/maxXXSize], White, 
    Point[{(maxXXSize - 1)/2, (maxYYSize - 1)/2}], PointSize[4/maxXXSize], 
    Black, Point[{(maxXXSize - 1)/2, (maxYYSize - 1)/2}]
    }];

you can be a bit more efficient, though, and Compile part of this, (don't worry too much about the "Formal" variables, I'm just using them in case you accidentally defined x, y, b, or c somewhere)
hh1C =
  With[{code = hh1[{[FormalX], [FormalY]}, [FormalB], [FormalC]]},
   Compile[{{[FormalX], _Complex}, {[FormalY], _Complex}, {[FormalB], 
_Real}, {[FormalC], _Real}}, 
    code,
    CompilationTarget -> "C"
    ]
   ];
newfunC[{{x_, y_}, n_}, b_, c_] := {hh1C[x, y, b, c], n + 1};

the speed up is relatively minor, but will increase with image size
Block[{
   maxIter = 63,
   maxXXSize = 500,
   maxYYSize = 500
   },
  unstablepic[0.125, -1.24];
  ] // AbsoluteTiming

{20.7144, Null}

Block[{
    newfun = newfunC,
    maxIter = 63,
    maxXXSize = 500,
    maxYYSize = 500
    },
   pic = unstablepic[0.125, -1.24];
   ] // AbsoluteTiming // First

18.9731

pic

Significant slowdown in a generalised implementation

One Answer

Add your own answers!

Ask a Question