7
$\begingroup$

Consider the following compiled function, which takes a $12 \times 5$ array $x_{ij}$ of real numbers and computes the triple sum $$ \sum_{k=1}^5 \sum_{i=1}^{12} \sum_{j=i+1}^{12} x_{ik} x_{jk}. $$

Compile[{{x, _Real, 2}},
    Sum[x[[i, k]]*x[[j, k]],
        {k, 5}, {i, 12}, {j, i + 1, 12}
    ]
]

When I attempt to evaluate this, Mathematica 10.1 on Windows throws the following error at me:

Compile::cpintlt: i+1 at position 2 of x[[i+1,5]] should be either a nonzero integer or a vector of nonzero integers; evaluation will use the uncompiled function. >>

It looks like the indices of the Sum[] are being inserted before the expression is passed to Compile[]. Why does this happen? Note that if I re-index the sum as $$ \sum_{k=1}^5 \sum_{i=1}^{12} \sum_{j=1}^{i-1} x_{ik} x_{jk} $$ Mathematica has no problem compiling the resulting function:

Compile[{{x, _Real, 2}},
    Sum[x[[i, k]]*x[[j, k]],
        {k, 5}, {i, 12}, {j, i-1}
    ]
]

This compiles just fine. Why does this version work, and not the previous version?

Improve this question
$\endgroup$
  • 3
    $\begingroup$ The original one also works if you convert the threefold sum into three nested sums. I believe you've encountered a bug, or at the least, a limitation of Compile (of which there are many). $\endgroup$ – Oleksandr R. Sep 28 '15 at 23:45
  • $\begingroup$ @m_goldberg That isn't true. Mathematica's default behavior in the case of an empty sum (lower index > upper index) is to return 0, even inside a CompiledFunction. You can try, for example Sum[x[[i]], {i, 13, 12}]. $\endgroup$ – David Zhang Sep 29 '15 at 1:20
  • $\begingroup$ I don't understand why you'd want to write the sums which contain combinations of indices that don't exist. The first sum should have $i$ going from 1 to 11 rather than 1 to 12 and the second sum should have $i$ going from 2 to 12. In addition for this example the sum can be greatly made more efficient (in terms of number of operations) to (Sum[(Sum[x[[i, k]], {i, 12}])^2, {k, 5}] - Sum[x[[i, k]]^2, {k, 5}, {i, 12}])/2 or (Total[Total[x]^2 - Total[x^2]])/2. This follows from expanding the square of the column sums. $\endgroup$ – JimB Sep 29 '15 at 3:43
  • $\begingroup$ @JimBaldwin You're absolutely right, and if I were actually using this code I would definitely adopt one of those simplifications. I used this only as an intentionally simplistic example to demonstrate a potentially buggy behavior I found. $\endgroup$ – David Zhang Sep 29 '15 at 4:10
  • 8
    $\begingroup$ The third argument of Compile might be what you need — e.g. Compile[{{x, _Real, 2}}, Sum[x[[i, k]]*x[[j, k]], {k, 5}, {i, 12}, {j, i + 1, 12}], {{i, _Integer}, {j, _Integer}, {k, _Integer}}]. $\endgroup$ – Stephen Luttrell Sep 29 '15 at 8:08
6
$\begingroup$

This seems rather be a error of determining the type of i, j and k than Sum itself. When you introduce your iterator variables in a Module and make clear that they are of type integer, it compiles fine for me:

fc = Compile[{{x, _Real, 2}},
  Module[{i = 0, j = 0, k = 0}, 
   Sum[x[[i, k]]*x[[j, k]], {k, 5}, {i, 12}, {j, i + 1, 12}]]
]

Btw, I want to note that StephenLuttrell's solution of giving the type as last argument works too.

Improve this answer
$\endgroup$

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.