7
$\begingroup$

I want to compile an expressions that contains nested Sum-expressions. A simple example that gives me problems is:

sumsum = Compile[{{n, _Integer}}, Sum[1/Sum[i + j, {j, n}], {i, n}]]

Note: This is not the problem I want to solve, it is just a simple example code, where the problem, which is described below, occurs.

While the expression compiles, using CompilePrint shows that the whole evaluation does not use compiled code, because there clearly is a call to MainEvaluate when the first call to Sum appears:

Needs["CompiledFunctionTools`"]
CompilePrint[sumsum]

A quick fix would be to use Plus@@Table instead of Sum, but this combination is in general slower than using Sum. Indeed

Plus @@ Table[1/Sum[i + j, {j, n}], {i, n}]

Does not show any call to MainEvaluate and runs faster than the uncompiled code.

Any suggestions on how to get nested sums to compile correctly? (once again, I am interested in a general solution to compiling nested sums, not a specific solution to the example stated above)

$\endgroup$
  • 1
    $\begingroup$ If this is your real problem, you can easily get an analytic solution. $\endgroup$ – b.gates.you.know.what Sep 24 '14 at 9:43
  • $\begingroup$ @b.gatessucks: This is of course not my real problem, but the real problem has the same structure. I edited the example slightly to show that there is a dependency on external data which leaves no analytic solution. $\endgroup$ – Wizard Sep 24 '14 at 9:48
  • 1
    $\begingroup$ Then a quick analysis shows that your nested Sum can be written as sumsum = Compile[{{data, _Real, 1}, {n, _Integer}}, Module[{t = Total@data[[;; n]]}, Total[(1/(n data[[;; n]] + t))]]] $\endgroup$ – xzczd Sep 24 '14 at 10:01
  • $\begingroup$ @xzczd: Thanks for the comment, but the problem above is not the problem I want to solve, it's just an example. The introduction of data seems to cause even more confusion, so I will edit it out again. Basically what I am trying to ask is how to compile nested Sums efficiently. $\endgroup$ – Wizard Sep 24 '14 at 11:17
  • $\begingroup$ Well, I believe the nested Sum can be generally avoided. BTW, according to my quick test, Plus@@Table is faster than Sum, and Total@Table will be even faster. $\endgroup$ – xzczd Sep 24 '14 at 11:29
12
$\begingroup$

Let me elaborate my comment into an answer. To make the nested Sum compiled, let's first have a close look at the compiling result of code containing one Sum:

sum = Compile[{{n, _Integer}}, 1/Sum[3.141 + j, {j, n}]];
Needs["CompiledFunctionTools`"]
p1 = CompilePrint@sum

enter image description here

It's not hard to notice that Sum is actually translated into a loop by Compile. It's reasonable to guess, if we express the summation with a loop instead of Sum, it will be finally translated into the same result, and it's true:

sum2 = Compile[{{n, _Integer}}, Module[{k = 0.}, 1/(Do[k += 3.141 + j, {j, n}]; k)]];
p2 = CompilePrint@sum2;
Grid[{{p1, p2}}, Frame -> All]

enter image description here

You can see the compiling results are almost same.

Manually replacing Sum with Module[…, Do[…]] is tedious, let's make use of code generation technique to facilitate the replacement:

sum[expr_, iter_] := Module[{k = 0.}, Do[k += expr, iter]; k];

sum2test = Hold@Compile[{{n, _Integer}}, 1/sum[3.141 + j, {j, n}]] //. DownValues@sum // 
   ReleaseHold;

sum2test@1000 == sum2@1000
(* True *)

Then since the Module[… Do[… ]] algorithm is closer to the the one adopted inside Compile, it's reasonable to guess that Compile can translate it more easily i.e. it may be compiled successfully in a situation that Sum fails, for example, nested Sum, and again, it's true:

sumsum2 = Hold@Compile[{{n, _Integer}}, sum[1/Sum[i + j, {j, n}], {i, n}]] //. 
DownValues@sum // ReleaseHold;

CompilePrint@sumsum2    

enter image description here

Pictured by Simon Wood's shadow.

Now the code is completely compiled and I believe it had the same performance as the compiled nested Sum (if it was correctly compiled).

Another acceptable solution is

sumsum3 = Compile[{{n, _Integer}}, Total@Table[1/Sum[i + j, {j, n}], {i, n}]];

It's a little faster than sumsum2:

sumsum3[10^4] // AbsoluteTiming
sumsum2[10^4] // AbsoluteTiming
{4.073000, 0.000109848}
{4.156000, 0.000109848}

However, what I really want to say is, I believe a nested Sum is not the fastest solution for your original problem at all, considering the facts that most of the arithmetic functions in Mathematica are Listable and taking advantage of the vectorization is one of the most important approaches to speed up codes even inside Compile. Here's a vectorized version of your toy code:

sumsum0 = Compile[{{n, _Integer}}, Total[1/(n Range@n + Total@Range@n)]];
sumsum0[10^4] // AbsoluteTiming
{0.001000, 0.000109848}
$\endgroup$
  • $\begingroup$ +1 - a very well-thought out and clear answer $\endgroup$ – dr.blochwave Sep 25 '14 at 7:26

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.