Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I have a function which calls other functions and I would like to compile it to C code. For example:

foo[p_] := 0.7 Exp[-p^2]
f[a_, n_] := 
 Module[{s, p, i}, s = 0; p = a; 
  For[i = 0, i < n, i++, s += foo[p]; p *= a;]; s]

Now if I do

cf = Compile[{{a, _Real}, {n, _Integer}}, f[a, n]]

and print the C code with

Needs["CompiledFunctionTools`"]
CompilePrint[cf]

the output is

R1 = MainEvaluate[ Hold[f][ R0, I0]]

which is obviously not desired.

When I try Compile[{{a, _Real}, {n, _Integer}}, Evaluate@f[a, n]], it just gives a function which returns 0...

share|improve this question
1  
Have a look at my answer here mathematica.stackexchange.com/q/24595/66 , it does what you want with the function CompileExpand and if the names of the user functions called are defined with upper cases. –  faysou Mar 19 at 7:42

1 Answer 1

up vote 12 down vote accepted

You need to inject the body of your function definition(s) into the Compile call without evaluating it. One possibility is (I had to change your Module a bit to correct a type error)

foo[p_] := 0.7 Exp[-p^2]
f[a_, n_] := Module[{s = 0.0, p = a, i},
  For[i = 0, i < n, i++, s += foo[p]; p *= a;]; s]

cf = ReleaseHold[Hold[Compile[{{a, _Real}, {n, _Integer}}, f[a, n]]] 
       /. DownValues[f] /. DownValues[foo]]

More explanation

First of all you need to understand that when you make a function definintion, Mathematica stores a replacement rule. There are different kinds of such rules and the most important ones are DownValues and OwnValues. OwnValues are used when you make assignments like x=1 or x:=1. DownValues are for function definitions. Let take a look at the rules of foo

DownValues[foo]
(* {HoldPattern[foo[p_]] :> 0.7 Exp[-p^2]} *)

Very sloppy speaking, when you make any calculation in Mathematica then it looks at your code and checks what rules are connected with a symbol. Then it uses these rules and makes the replacement, if the pattern matches.

In your case, the problem with your Compile call is that

  1. Compile has the attribute HoldAll which prevents that f[a,n] is replaced with its definition.
  2. If you force the evaluation of the expression f[a,n] using Evaluate then it does too much, because it doesn't only insert the definition body of f (the Module part) but it evaluates it right away. This is not good because at this point, you don't have numeric values for a and n and therefore, the For loop doesn't do anything at all.

At this point, it is clear what you want: You have made a definition for f[a,n], how can you get the fully expanded function body? Two ingredients: First, you use the definition rules (DownValues) by yourself and second, you prevent the evaluation using Hold. There you go..

Hold[f[a, n]] /. DownValues[f]
(* Hold[Module[{s$ = 0., p$ = a, i$}, 
         For[i$ = 0, i$ < n, i$++, s$ += foo[p$]; p$ *= a;]; s$]] 
*)

As you see we are not finished, because the definition of foo was not expanded, but this is not a problem, because we can do it the same way. All this combined with Compile gives the above solution.

share|improve this answer
    
Would you care to explain a bit how this works? What is the purpose of DownValues? –  Danvil Mar 19 at 9:15
    
@Danvil I have added some information in my answer. –  halirutan Mar 19 at 10:19
    
Thank you very much! So am I assuming correctly that Mathematica can not compile a recursive function to C? E.g. f[n_] := If[n == 0, 1, n f[n - 1]]. –  Danvil Mar 19 at 13:52
    
@Danvil Yes, I discussed this issue several times with e.g. @Oleksandr R. who has quite some knowledge about such stuff. When I remember correctly, then there is no general way to make a recursive compiled function. This goes along with Leonids answer here. –  halirutan Mar 19 at 21:00

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.