3
$\begingroup$

I have a piece of code that generates a figure, performing a calculation to generate each pixel. At a certain resolution it runs fine, but after a high enough resolution it fails. Mathematica beeps (sounds more line a 'dong') and just stops after a few seconds with no error warnings. How can I run some diagnostics on my code to try and figure out the problem?

I assume it's a memory problem, but when I use the MaxMemoryUsed[] command it was only ~20Mb which doesn't seem that much.

The code is shown below. I adjust the parameter 'res'. Up to res=12 it runs fine but thereafter it fails.

Clear["Global`*"];
startTime = AbsoluteTime[];
f[{rx_, ry_}] := 
 Module[{e1, e2, e3, d1, d2, d3, a, b, c, sol1, sol2, data, round, 
   rArray},
  {e1, e2, e3} = # & /@ {{0, -1}, {Sqrt[3]/2, 1/2}, {-Sqrt[3]/2, 1/2}};
  If[
   Norm[{rx, ry} - e1] <= 1/3 || Norm[{rx, ry} - e2] <= 1/3 || 
    Norm[{rx, ry} - e3] <= 1/3,
   Null,
   {d1, d2, d3} = # - {rx, ry} & /@ {e1, e2, e3};
   {a, b, c} = Norm[#]^-3 & /@ {d1, d2, d3};
   round = 0.000001;
   sol1 =
    If[0 <= ((b + c)^2 - a^2)/(4 b c) <= 1, 0,
     Re[{x, y} /. Solve[
        Round[a b (d2 - d1), round] Sin[{x, y}.(d1 - d2)]
          + Round[b c (d3 - d2), round] Sin[{x, y}.(d2 - d3)]
          + Round[a c (d3 - d1), round] Sin[{x, y}.(d1 - d3)]
         == 0, {x, y}]]];
   sol2 = 
    Sqrt[a^2 + b^2 + c^2 + 2 a b Cos[#.(d1 - d2)] + 
        2 b c Cos[#.(d2 - d3)] + 2 a c Cos[#.(d1 - d3)]] & /@ sol1;
   Min[sol2]
   ]
  ]
res = 20;
rArray = Table[{rx, ry}, {rx, -Sqrt[3]/2, Sqrt[3]/2, 
    Sqrt[3]/res}, {ry, 3/4, -3/4, -3/(2 res)}];
diagram = Map[f, rArray, {2}];
ArrayPlot[diagram, ColorFunction -> "Rainbow", 
 ColorRules -> {0 -> White, Null -> Gray}]
endTime = AbsoluteTime[] - startTime
$\endgroup$
  • $\begingroup$ In Mma v10 this code runs without any problem. In Mma v9 the kernel crashes while executing diagram = Map[f, rArray, {2}]; $\endgroup$ – Karsten 7. Jul 23 '14 at 12:41
  • $\begingroup$ Does menu command Help > Why the Beep? tell you anything? $\endgroup$ – Mr.Wizard Jul 23 '14 at 12:44
  • $\begingroup$ @Mr.Wizard : "The kernel Local has quit (exited) during the course of an evaluation." $\endgroup$ – Tom Jul 23 '14 at 12:46
  • 1
    $\begingroup$ @Karsten7. : thanks for the info. Could be a good excuse to get my supervisor to upgrade! =P $\endgroup$ – Tom Jul 23 '14 at 12:46
  • $\begingroup$ What's the use of the Null in the If function? $\endgroup$ – Sjoerd C. de Vries Jul 23 '14 at 12:51
3
$\begingroup$

Get this code working

You can get the code running under Mma v9, by adding options to Solve:
Adding

VerifySolutions -> False

should work.
And also

VerifySolutions -> True, WorkingPrecision -> MachinePrecision

works on my PC, if res isn't too big.
For bigger res (tested for 50)

VerifySolutions -> False, WorkingPrecision -> MachinePrecision

will make the code working.


Parallelization

You didn't ask for this, but if you execute

ParallelEvaluate[Off[Solve::ifun]]

before the posted code, the Solve::ifun messages will be switched off.
And replacing Map with ParallelMap:

diagram=ParallelMap[f,rArray,{2}];

will speed up the code significantly, if you have more than one core.

$\endgroup$
  • 1
    $\begingroup$ Brilliant! This is working for me and also you've cut my computation time in half, which is great because I'm doing a lot of this sort of thing and often I leave things computing for hours, so thanks! $\endgroup$ – Tom Jul 26 '14 at 9:08

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.