Skip to main content

compile - Catching and debugging numerical errors in compiled functions


I have a compiled function (that was created using the second method of this question). Basically it's using a compiled closure. The problem is that for some parameter values it does not evaluate. Here is a small example:


(*just some function that has a singularity at x=5*)
f = Compile[{{x, _Real, 0}, {a, _Real, 0}}, Sin[a*x*Pi/180]/(x - 5)
, RuntimeOptions -> {"Speed", "EvaluateSymbolically" -> False}];


(*The wrapper just increments the x values by diff. Aa is the
variable that is supposed to be inlined in the closure*)
Block[{Aa},
With[{wrapper =
Compile[{{x, _Real, 0}, {diff, _Real, 0}},
Evaluate@f[x + diff, Aa]]},
min = Compile[{{xs, _Real, 1}, {Aa, _Real, 0}, {diff, _Real, 0}},
First@Sort[wrapper[#, diff] & /@ xs]
, CompilationOptions -> {"InlineCompiledFunctions" -> True,

"InlineExternalDefinitions" -> True,
"ExpressionOptimization" -> True},
RuntimeOptions -> {"Speed"}](*compile*)
]
];

(*The warning is harmeless*)
(*CompiledFunction::cfsa:Argument diff+x at position 1 should be a machine-size real number. >>*)

(*The function is compiled*)

Needs["CompiledFunctionTools`"]
CompilePrint[min]


3 arguments
9 Integer registers
11 Real registers
3 Tensor registers
Underflow checking off
Overflow checking off

Integer overflow checking off
RuntimeAttributes -> {}

T(R1)0 = A1
R0 = A2
R1 = A3
I7 = 0
I6 = -5
I3 = 1
R6 = 3.141592653589793

I5 = 180
Result = R7

1 I2 = Length[ T(R1)0]
2 I8 = I7
3 T(R1)1 = Table[ I2]
4 I4 = I7
5 goto 20
6 R3 = GetElement[ T(R1)0, I4]
7 R2 = R1

8 R5 = R2 + R3
9 R4 = R0
10 R8 = I5
11 R9 = Reciprocal[ R8]
12 R8 = R6 * R9
13 R9 = R4 * R5 * R8
14 R8 = Sin[ R9]
15 R9 = I6
16 R10 = R5 + R9
17 R9 = Reciprocal[ R10]

18 R8 = R8 * R9
19 Element[ T(R1)1, I8] = R8
20 if[ ++ I4 < I2] goto 6
21 T(R1)2 = Sort[ T(R1)1]]
22 R7 = Part[ T(R1)2, I3]
23 Return

But then here are the use cases:


(*everything works OK even though the singularity is in the arguments*)
min[{1,5,5},1,0]

(*=> -0.0043631*)

But for some arguments it does not:


(*But if all the arguments are singularites: BLAM!*)
min[{5,5,5},1,0]

(*=> CompiledFunction[{x,a},Sin[(a x \[Pi])/180]/(x-5),-CompiledCode-][5,Aa] *)


   CompiledFunction::cfne: Numerical error encountered; proceeding with uncompiled evaluation. >>

CompiledFunction::cfse: Compiled expression Aa should be a machine-size real number. >>
CompiledFunction::cfex: Could not complete external evaluation at instruction 2; proceeding with uncompiled evaluation. >>

Questions:



  • Why is this working if some of the arguments give division by zero, but does not evaluate if all of the arguments are incorrect?

  • How can one debug such problems? There is nothing obviously wrong at instruction 2 in the compiled output...

  • How to guard against and/or catch such errors in compiled code?




Comments

Popular posts from this blog

plotting - How to draw lines between specified dots on ListPlot?

I would like to create a plot where I have unconnected dots and some connected. So far, I have figured out how to draw the dots. My code is the following: ListPlot[{{1, 1}, {2, 2}, {3, 3}, {4, 4}, {1, 4}, {2, 5}, {3, 6}, {4, 7}, {1, 7}, {2, 8}, {3, 9}, {4, 10}, {1, 10}, {2, 11}, {3, 12}, {4,13}, {2.5, 7}}, Ticks -> {{1, 2, 3, 4}, None}, AxesStyle -> Thin, TicksStyle -> Directive[Black, Bold, 12], Mesh -> Full] I have thought using ListLinePlot command, but I don't know how to specify to the command to draw only selected lines between the dots. Do have any suggestions/hints on how to do that? Thank you. Answer One possibility would be to use Epilog with Line : ListPlot[ {{1, 1}, {2, 2}, {3, 3}, {4, 4}, {1, 4}, {2, 5}, {3, 6}, {4, 7}, {1, 7}, {2, 8}, {3, 9}, {4, 10}, {1, 10}, {2, 11}, {3, 12}, {4, 13}, {2.5, 7}}, Ticks -> {{1, 2, 3, 4}, None}, AxesStyle -> Thin, TicksStyle -> Directive[Black, Bold, 12], Mesh -> Full, Epilog -> { Line[ ...

list manipulation - Selecting multiple columns from a matrix?

Sample data: data = { {{2013, 1, 1}, 24.13, 167.67, 231.82}, {{2013, 1, 2}, 32.15, 170.92, 225.99}, {{2013, 1, 3}, 35.43, 172.68, 221.67}, {{2013, 1, 4}, 36.73, 173.05, 218.32}, {{2013, 1, 5}, 58.19, 165.96, 197.05}, {{2013, 1, 6}, 69.99, 163.50, 187.52}, {{2013, 1, 7}, 71.37, 154.21, 175.58}, {{2013, 1, 8}, 72.51, 149.66, 163.25}}; I want a DateListPlot with three graphs, so for a matrix formed by columns 1 and 2, one for columns 1 and 3, and 1 for columns 1 and 4. At the moment I'm using this code: data2 = Transpose[{data[[All, 1]], data[[All, 2]]}]; data3 = Transpose[{data[[All, 1]], data[[All, 3]]}]; data4 = Transpose[{data[[All, 1]], data[[All, 4]]}]; DateListPlot[{data2, data3, data4}, Joined -> True, Filling -> {3 -> {1}}] but I have a hunch that this can be done more efficiently. I don't like the Transpose s in particular. Any ideas? edit (for extra credit) What if I need to multiply the second column by 2, which in my solution is simp...

plotting - How to draw some 2D curves in a 3D plot?

I want to plot several 2D curves ( Sin[x] at y=0 , Cos[x] at y=1 ) in a 3D plot, what function should i use? f[x_, y_] := Piecewise[{{0, y 0}, {Sin[x], y == 0}}] g[x_, y_] := Piecewise[{{0, y 1}, {Cos[x], y == 1}}] Plot3D[{f[x, y], g[x, y]}, {x, 0, 10}, {y, 0, 10}] doesn't work. Any ideas? Answer ParametricPlot3D[{{Sin[u], 0, u}, {Cos[u], 1, u}}, {u, 0, 20}, BoxRatios -> 1, PlotRange -> {{-2, 2}, {-.5, 2}, {0, 20}}] === update - general thoughts === I've seen the comments to this question and a few other approaches to the waterfall (or wire) plot are given in this question: Plotting several functions . Also I cannot not mention one of my favorite Mathematica blogs - Archery with this post on the story of discovery of 1st pulsar and an its infamous graphics - code is given in the post: