error when defining a compiled function with ReplacePart

Question

I always have difficulty in writing a non-trivial compiled function. I'm using Mathematica 9. Please see the following code

Clear[longtermCompiled]; 
longtermCompiled = 
 Compile[{{M, _Real, 2}}, 
  Module[{len = Length[M]}, 
   ReplacePart[Transpose[M] - IdentityMatrix[len], {len, i_} -> Table[1, {len}]]]];

where M is expected to be a matrix (two-dimensional). Upon executing the above definition, it complained:

ReplacePart::argrx: ReplacePart called with 2 arguments; 3 arguments are expected.

I don't understand this error message, because the following code works without any problems:

Clear[longtermX]; 
longtermX[M_] := 
 Module[{len = Length[M]}, 
  ReplacePart[Transpose[M] - IdentityMatrix[len], {len, i_} -> Table[1, {len}]]];

What's wrong?

Answer 1

It's not you, it's Mathematica. You are not expected to know this, but basically, in compiled code, ReplacePart merely acts as syntactic sugar for setting a part, i.e.:

l = Range[3]; ReplacePart[l, 2 -> 0]
(* -> {1, 0, 3} *)

would be compiled (but see below) into exactly the same bytecode as

l = Range[3]; Block[{l = l}, l[[2]] = 0; l]
(* -> {1, 0, 3} *)

(which you can easily verify if you so desire).

It should be obvious from this equivalence that the pattern syntax allowed by ReplacePart in top-level code cannot be used in compiled code; indeed, no patterns or even Rules are accepted by the compiler. (As such, the example given above actually will not compile as written.) To get around the latter limitation it is possible to use an undocumented (since v6, but see version 5.2) variant syntax of ReplacePart which works both at the top level and in compiled code: replace the Rule with a Sequence, and reverse the order of the arguments. That is, the first example now looks like this:

l = Range[3]; ReplacePart[l, 0, 2]
(* -> {1, 0, 3} *)

So, this should suffice to explain (a) that ReplacePart can only be used in a strictly limited sense within compiled code and (b) the meaning of the strange message that you receive as a result of your attempt to use it in an unsupported way.

It happens in this case that you can re-write your ReplacePart call in terms of Part and Set:

Clear[longtermX];
longtermX[M_?MatrixQ] := Module[{len = Length[M], tmp},
   tmp = Transpose[M] - IdentityMatrix[len];
   tmp[[len, All]] = Sequence@ConstantArray[1, len];
   tmp
  ];

However, this does not really help in practice because a ragged array (which the VM does not support) is produced as output, and this usage of Sequence together with Set/Part is not properly understood by the compiler anyway. What one might try instead, and which will work in compiled code, is something along the lines of tmp[[len, All]] = 1 in place of tmp[[len, All]] = Sequence@ConstantArray[1, len], then tmp2 = ConstantArray[Last[tmp], len]; tmp = Most[tmp] to have the two (full-rank) pieces of the array separately.