Error from Compile with a sum over two indices

Question

I have the following code:

compVDiag = Compile[
  {{Nn, _Integer}, {M, _Integer}, {n, _Integer}, {basisvecs, _Integer,
     2}, {c1m, _Integer, 2}, {c2m, _Integer, 2}},
  Block[
   {zpow, wpow, zlen, wlen},
   zpow = Union[basisvecs[[n, 1 ;; Nn]]];
   wpow = Union[basisvecs[[n, Nn + 1 ;; Nn + M]]];
   zlen = Length[zpow];
   wlen = Length[wpow];

      Sum[
        2^(-zpow[[i]] - zpow[[j]])*c1m[[n, zpow[[i]] + 1]]*
         c1m[[n, zpow[[j]] + 1]]*
         Binomial[zpow[[i]] + zpow[[j]], zpow[[i]]], {i, 1, 
         zlen - 1}, {j, i + 1, zlen}]
       +
       Sum[
        2^(-wpow[[i]] - wpow[[j]])*c2m[[n, wpow[[i]] + 1]]*
         c2m[[n, wpow[[j]] + 1]]*
         Binomial[wpow[[i]] + wpow[[j]], wpow[[i]]], {i, 1, 
         wlen - 1}, {j, i + 1, wlen}]

  ]

which gives me some warnings:

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

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

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

General::stop: Further output of Compile::cpintlt will be suppressed during this calculation. >>

CompilePrint shows the following weird snippet:

11  R2I4T(I1)3I8T(I1)4 = MainEvaluate[                                 \
                                                                      \
                                                                      \
                                    -zpow[[1]] - zpow[[i + 1]]
Function[{Nn, M, n, basisvecs, c1m, c2m, zlenCompile$1, \
    wpowCompile$2, wlenCompile$3, zpowCompile$4}, Block[{zlen = \
zlenCompile$1, wpow = wpowCompile$2, wlen = wlenCompile$3, zpow = \
    zpowCompile$4}, {Length[{2                           c1m[[n,zpow[[1]] \
+ 1]] c1m[[n,zpow[[i + 1]] + 1]] Binomial[zpow[[1]] + zpow[[i + 1]], \
zpow[[1]]]}[[1]]], zlen, wpow, wlen, zpow}]][ I0, I1, I2, T(I2)0, \
T(I2)1, T(I2)2, I4, T(I1)3, I8, T(I1)4]]

Could I get some help in understanding this? Seems like i is not assigned a value yet when a term in the sum is calculated, even though j starts at i+1? How can I remedy this?

Answer 1

It happens quite often that compilation of a correctly running Mathematica function for some or other reason fails. The reason for that is not always easy to trace. After some testing I feel that the problem you ran into is close to a bug.

First let us further simplify your problem, showing the same warnings:

fc=Compile[{}, Block[{mat},
  mat=Array[1&, {10,10}];
  Sum[mat[[i,j]], {i,1,10},{j,i,10}]]];

(*
During evaluation of In[75]:= Compile::cpintlt: i at position 3 of mat[[1,i]] should be either a nonzero integer or a vector of nonzero integers; evaluation will use the uncompiled function. >>
During evaluation of In[75]:= Compile::cpintlt: i at position 3 of mat[[1,i]] should be either a nonzero integer or a vector of nonzero integers; evaluation will use the uncompiled function. >>

*)

These two warnings are strange: there is no position 3 in mat[[1,i]], the variable i indeed is an integer, and there is no expression mat[[1, i]] in the code.

The warnings have nothing to do with the double summation. When we replace mat[[i,j]] by 1 or i+j, we get no warnings and there is no call of MainEvaluate.

The warnings only turn up when the second underlimit depends on the first index. For example,

fc=Compile[{}, Block[{mat},
  mat=Array[1&, {10,10}];
  Sum[mat[[i,j]], {i,1,10},{j,1,i}]]];

compiles nicely.

Somehow, your problem is a combination of using a Part-function in the first argument of your double sum and a dependency of the underlimit of the second index of the first index. In all other cases it works, and I think it should work in this case as well. So very likely a bug.

When we replace the double summation by a nested summation, there are no warnings:

fc=Compile[{}, Block[{mat},
  mat=Array[1&, {10,10}];
  Sum[Sum[mat[[i,j]], {j,i,10}], {i,1,10}]]];

But CompilePrint shows a call of the kernel for evaluating the nested sums, and moreover the inner sum is not recognized as returning an integer; the result is now a real:

fc[]
(* 55. *)

So this does not help.

When you want to use the double summation, the only workaound I see is to formulate it in such a way that the second underlimit does not depend on the first index:

fc=Compile[{}, Block[{mat},
  mat=Array[1&, {10,10}];
  Sum[mat[[i,j]], {i,1,10}, {j,10,i,-1}]]];

Ridiculous, but it works. CompilePrint does not show a MainEvaluate.

fc[]
(* 55 *)