6

(Disclaimer: I realize that the addition of matrices a and b below can be done very fast just by typing a+b. The code below is just meant to illustrate behavior similar to what I observe in much more complicated code, too lengthy to reproduce here).

Suppose I want to add matrices a and b

n = 50;

I can do it via a double loop in a single kernel, this takes 0.02 seconds:

a = RandomReal[{}, {n, n}];
b = RandomReal[{}, {n, n}];
c1 = RandomReal[{}, {n, n}];
Do[
  Do[
   c1[[i, j]] = b[[i, j]] + a[[i, j]],
   {j, 1, n}],
 {i, 1, n}] // AbsoluteTiming

Out[363]= {0.0240014, Null}

I can try to parallelize the outer loop -- I use 4 kernels -- this takes 25 seconds -- ouch -- presumably due to the fact that the master kernel is in constant communication with the various kernels, updating c2.

SetSharedVariable[c2];
c2 = RandomReal[{}, {n, n}];
ParallelDo[
  Do[
   c2[[i, j]] = b[[i, j]] + a[[i, j]],
   {j, 1, n}],
  {i, 1, n}] // AbsoluteTiming

Out[366]= {25.8984813, Null}

I can try to reduce communication to an absolute minimum, by introducing a temporary array on each kernel. Surprisingly, this still requires 0.3 seconds. 

SetSharedVariable[c3];
UnsetShared[vtemp];
c3 = RandomReal[{}, {n, n}];
ParallelDo[
  vtemp = Table[0., {n}];
  Do[
   vtemp[[j]] = b[[i, j]] + a[[i, j]],
   {j, 1, n}];
  c3[[i]] = vtemp,
  {i, 1, n}, Method -> "CoarsestGrained"] // AbsoluteTiming

Out[370]= {0.2940168, Null}

Again, my own code is far more complicated. But it shares with this example the need to distribute the task of updating a piece of a matrix. Updating the global variable/matrix invariably takes more time than executing the code on a single processor, even when the calculations on each processor are far more complex that in the above. Are there techniques to mitigate this problem?

Eli Lansey
  • 7,499
  • 3
  • 36
  • 73
Eric
  • 333
  • 1
  • 10
  • 1
    Is there any chance you can avoid using Set, or any other functions with side effects? Please see here: http://mathematica.stackexchange.com/a/1771/12 Since your problem is parallelizable, perhaps it can be formulated in a way that Mathematica can Parallelize automatically. Note that Parallelize will handle a lot more than Map, Do & Table. It can also deal with e.g. Inner, Outer, MapThread, etc. (E.g. this matrix addition is easily reformulated in a functional way using MapThread) – Szabolcs Feb 16 '12 at 16:50
  • Actually a quick experiment has shown that MapThread[Plus, {a, b}]; // AbsoluteTiming is much faster than Parallelize[MapThread[Plus, {a, b}]]; // AbsoluteTiming (regardless of the Method setting), so my advice from the previous comment is not very good ... – Szabolcs Feb 16 '12 at 16:57

1 Answers1

5

Is there any reason you can't directly construct the Table, as in:

DistributeDefinitions[a,b]
c1 = ParallelTable[b[[i, j]] + a[[i, j]], {j, 1, n}, {i, 1, n}]

On my machine, with n=500, Table[b[[i, j]] + a[[i, j]], {j, 1, n}, {i, 1, n}];// AbsoluteTiming takes 0.59 s, while ParallelTable takes only 0.16 s. With n=1000, Table takes 2.1s while ParallelTable takes only 0.47s.

For future reference, unlike in MATLAB or other languages, in Mathematica you do not need to first "initialize" the table, and then fill it.

Eli Lansey
  • 7,499
  • 3
  • 36
  • 73
  • It is not a bad idea to initialize large lists, even in Mathematica. See this post by Leonid, item 3.1 – rm -rf Feb 16 '12 at 17:27
  • @R.M I've seen that amazing post. I'd love to know the constraints on that, because IME using Set is always slower than using Table or Map (sometimes those methods are faster than Reap/Sow, too.) – Eli Lansey Feb 16 '12 at 17:30
  • OK, problem is that in my real code, c aleady exists and is being updated in the loop. A better example might have been c = c + a as opposed to c = a + b. This necessitates having c available to the kernels. – Eric Feb 16 '12 at 19:43
  • To make sure I understand, given your code above, you want something like c[[i,j]] = c[[i,j]] + b[[i,j]] in the Do loop? – Eli Lansey Feb 16 '12 at 20:25
  • And, if that's the case, any reason you can't do d[[i,j]] = c[[i,j]] + b[[i,j]], followed by c=d after updating d? – Eli Lansey Feb 16 '12 at 20:32