Using Apply inside Compile

14

7

In this great answer a compiled version of the Nelder-Mead algorithm is presented.

Since it works on arbitrary dimensions (i.e. arbitrary number of arguments), it has to use apply on the objective function. The problem is that Apply is not directly supported inside Compile. To overcome this limitation the following code is used:

(* Produces compiled code for the Nelder-Mead algorithm with the objective function inlined *)
ClearAll[apply];
SetAttributes[apply, HoldAll];
apply[f : (_Function | _CompiledFunction), vars : {__Symbol}] :=
    With[{applied := f @@ vars},
      Function[arglist, Block[vars, vars = arglist; applied]]
     ];

This seems to work (inside the package), but I have no idea, how it works.

Could somebody explain the techniques behind this code snippet? Specifically:

  • What does SetDelayed do inside With?
  • Why does vars appear twice in the arguments to Block?

  • Why doesn't this minimal example work? (I believe it mimics what is done in the mentioned answer...)

    Clear[a, b, c, x, y, z, objectiveFunction, cfunc];
    objectiveFunction = Compile[{a, b, c, x, y, z} ,
       (a - x)^2 + 50 (b - y)^2 + (c - z)^2];
    
    cfunc = With[{f = apply[objectiveFunction, {a, b, c, x, y, z}]},
      Compile[ {{pts, _Real, 1}}, f@pts]
      ]
    << CompiledFunctionTools`
    CompilePrint[cfunc]
    

EDIT:

Here is the working snippet. Note that "InlineCompiledFunctions" must be set to True to avoid the call to MainEvaluate.

ClearAll[a, b, c, x, y, z, apply, objectiveFunction, cfunc];
(* Produces inlinable code for use inside Compile (where Apply is not \
supported directly) *)
SetAttributes[apply, HoldRest];
apply[f : (_Function | _CompiledFunction), vars : {__Symbol}] :=
       Function[arglist, Block[vars, vars = arglist; f @@ vars]];

objectiveFunction = Compile[{a, b, c, x, y, z} ,
   (a - x)^2 + 50 (b - y)^2 + (c - z)^2];
apply[objectiveFunction, {a, b, c, x, y, z}]
cfunc = With[{f = apply[objectiveFunction, {a, b, c, x, y, z}]},
  Compile[ {{pts, _Real, 1}}, f@pts,
    CompilationOptions -> {"InlineCompiledFunctions" -> True}]
  ]
<< CompiledFunctionTools`
CompilePrint[cfunc]

Ajasja

Posted 2012-05-04T08:19:13.547

Reputation: 13 114

This makes me wonder why is the a system option called ApplyCompileLength and why is it set to Infinity – Szabolcs – 2012-05-04T09:12:19.860

2@Szabolcs This is an entirely different story. Apply as discussed here evaluates at (pre) Compile-time and is a tool of code generation. The option you mention governs the auto-compilation of Apply, when it would be explicitly present in the Compile body when auto-compiling. Since Apply is only supported for 3 heads (List, Plus and Times) by Compile, it does not generally make sense to set the mentioned option to a finite number. In case you have one of those heads, you can of course reset it. – Leonid Shifrin – 2012-05-04T09:18:32.393

@Leonid (It was really a question, not a comment...) Is it documented somewhere that it's supported for these three heads? Do you know why that option is set to Infinity? – Szabolcs – 2012-05-04T09:27:09.497

@Szabolcs It is documented, but I can't tell you where off hand. For the option setting, Infinity is equivalent to saying "don't attempt to auto-compile at all". The fact that only 3 heads are supported means that Compile does not generally support Apply, which is probably the reason for this setting. For those heads, it can be reset manually, as I mentioned. – Leonid Shifrin – 2012-05-04T09:33:56.067

Answers

12

First question : With accepts a syntax like

  With[{var:=value}, expression]

in which case, value is injected into expression unevaluated. As far as I know, this syntax is not documented. You can achieve a similar effect with the replacement rules, by using

Unevaluated[expression]/.HoldPattern[var]:>value

There are some subtle differences between the semantics of With and repalcement rules though, mostly related to the treatment of nested scoping constructs and variable name conflicts in them.

Second question: vars appear twice because they must first be Block-ed, and then there is a massive assignment to them performed in the body of the Block. This is probably the most economical way of blocking a number of variables and assigning values to them simultaneously - otherwise a more complex code-generation will be needed. You can see another example of that in this answer (and if you look at the revision history for that answer, you can find an alternative, harder way to do this, in one of the previous revisions).

Third question: this does not work because the apply function was made HoldAll (which isn't quite necessary), and the pattern-matching does not work. There were some past discussions on this topic on SO, but can't find them right now. But I discussed this topic at length also in my book. The idea is that at the pattern-matching time, all seen by apply is a variable objectiveFunction, and because it does not evaluate it, the pattern _CompiledFunction is not matched. The solution is to make apply HoldRest, and then it works:

ClearAll[apply];
SetAttributes[apply, HoldRest];
apply[...]:=...

Leonid Shifrin

Posted 2012-05-04T08:19:13.547

Reputation: 108 027

Aha, so this works by wrapping the original function in another function which takes just a single argument (the list) and then mass assigning back to individual arguments? When is the Apply used in With evaluated? – Ajasja – 2012-05-04T08:55:00.537

Or to put it another way, why is With even necessary, since apply[f : (_Function | _CompiledFunction), vars : {__Symbol}] := Function[arglist, Block[vars, vars = arglist; f @@ vars]]; seems to work as well. – Ajasja – 2012-05-04T08:57:38.163

@Ajasja Only inside Block, and that's the entire point, because by that time, Block guarantees that the variables have no global values. – Leonid Shifrin – 2012-05-04T08:58:43.310

@Ajasja With isn't necessary. Your function is exactly equivalent. I suppose it was intended to make the code more readable. – Leonid Shifrin – 2012-05-04T09:01:44.023

Thanks, this was very instructive. – Ajasja – 2012-05-04T09:04:17.203

@Ajasja Glad to help, and thanks for the accept. – Leonid Shifrin – 2012-05-04T09:12:11.933

Leonid, many thanks for answering this, especially as it was arguably my responsibility to do so but I don't have time at the moment. I owe you one! (BTW, the With is just a remnant of a previous iteration of this function. Since it caused confusion, I should probably remove it in the posted answer.) – Oleksandr R. – 2012-05-04T10:21:12.883

3@OleksandrR. No problem :) You use some really nice stuff, which is a lot of fun to read and discuss. By the way, I take this opportunity to say that I am (among many others) quite impressed by that answer of yours which originated this question. I think we need more answers like that one, to show that one can often use Mathematica as an implementation language for some missing functionality, and that the absence or limitations of some built-in function don't mean that "Mathematica can not do this", or "can not do this fast enough". I think you did a great job proving such statements wrong. – Leonid Shifrin – 2012-05-04T10:44:31.980

1That SetDelayedish syntax in With is really useful. It's a shame it's not documented. +1! – Pillsy – 2012-05-04T18:46:41.160

@Pillsy Yes, it is useful, I agree. – Leonid Shifrin – 2012-05-04T18:47:53.167