Skip to main content

performance tuning - Benchmarking expressions


It often happens that I'd like to know which of two or more alternative expressions (that should evaluate to the same value) is fastest. For example, below are three ways to define a function:



fooVersion1[x_]/; x > 0 := E^(-1/x);
fooVersion1[x_] := 0;

fooVersion2[x_] := If[x > 0, E^(-1/x), 0];

fooVersion3 = Function[x, If[x > 0, E^(-1/x), 0]];

(Granted, strictly speaking, these three versions are not completely equivalent, but at least, for all real scalars x, fooVersion1[x], fooVersion2[x], and fooVersion3[x] should evaluate to the same value.)


Is there a way to determine which one is fastest that is both simple and reasonably reliable?



Answer




New method


Mathematica 10 introduced a Benchmark and BenchmarkPlot functions in the included package GeneralUtilities. The latter automates and extends the process described in the Legacy section below. Version 10.1 also introduced RepeatedTiming which is like a refined version of timeAvg that uses TrimmedMean for more accurate results. Here is an overview of these functions.


Before using BenchmarkPlot one must load the package with Needs or Get. Its usage message reads:


enter image description here


Benchmark has a similar syntax but produces a table of numeric results rather than a plot.


The Options and default values for BenchmarkPlot are:


{
 TimeConstraint -> 1.`,
 MaxIterations -> 1024, 
 "IncludeFits" -> False,

 "Models" -> Automatic
}

Notes:




  • BenchmarkPlot uses caching therefore if a function is changed after benchmarking its timings may not be accurate in future benchmarking. To clear the cache one may use:


    Clear[GeneralUtilities`Benchmarking`PackagePrivate`$TimingCaches]



  • There is a bug in 10.1.0 regarding BenchmarkPlot; a workaround is available.




Legacy method


I described my basic method of comparative benchmarking in this answer:
How to use "Drop" function to drop matrix' rows and columns in an arbitrary way?


I shall again summarize it here. I make use of my own modification of Timo's timeAvg function which I first saw here. It is:


SetAttributes[timeAvg, HoldFirst]

timeAvg[func_] := Do[If[# > 0.3, Return[#/5^i]] & @@ Timing@Do[func, {5^i}], {i, 0, 15}]


This performs enough repetitions of a given operation to exceed the threshold (0.3 seconds) with the aim getting a sufficiently stable and precise timing.


I then use this function on input of several different forms, to try to get a first order mapping of the behavior of the function over its domain. If the function accepts vectors or tensors I will try it with both packed arrays and unpacked (often unpackable) data as there can be great differences between methods depending on the internal optimizations that are triggered.


I build a matrix of timings for each data type and plot the results:


funcs = {fooVersion1, fooVersion2, fooVersion3};

dat = 20^Range[0, 30];

timings = Table[timeAvg[f @ d], {f, funcs}, {d, dat}];


ListLinePlot[timings]

Mathematica graphics


Here with Real data as a second type (note the N@):


timings2 = Table[timeAvg[f @ d], {f, funcs}, {d, N@dat}];

ListLinePlot[timings2]

Mathematica graphics


One can see there is little relation between the size of the input number and the speed of the operation with either Integer or machine-precision Real data.



The sawtooth shape of the plot lines is likely the result of insufficient precision (quantization). It could be reduced by increasing the threshold in the timeAvg function, or better in this case to time multiple applications of each function under test in a single timeAvg operation.


Comments

Post a Comment

Popular posts from this blog

functions - Get leading series expansion term?

Given a function f[x] , I would like to have a function leadingSeries that returns just the leading term in the series around x=0 . For example: leadingSeries[(1/x + 2)/(4 + 1/x^2 + x)] x and leadingSeries[(1/x + 2 + (1 - 1/x^3)/4)/(4 + x)] -(1/(16 x^3)) Is there such a function in Mathematica? Or maybe one can implement it efficiently? EDIT I finally went with the following implementation, based on Carl Woll 's answer: lds[ex_,x_]:=( (ex/.x->(x+O[x]^2))/.SeriesData[U_,Z_,L_List,Mi_,Ma_,De_]:>SeriesData[U,Z,{L[[1]]},Mi,Mi+1,De]//Quiet//Normal) The advantage is, that this one also properly works with functions whose leading term is a constant: lds[Exp[x],x] 1 Answer Update 1 Updated to eliminate SeriesData and to not return additional terms Perhaps you could use: leadingSeries[expr_, x_] := Normal[expr /. x->(x+O[x]^2) /. a_List :> Take[a, 1]] Then for your examples: leadingSeries[(1/x + 2)/(4 + 1/x^2 + x), x] leadingSeries[Exp[x], x] leadingSeries[(1/x + 2 + (1 - 1/x...
Post a Comment
Read more

mathematical optimization - Minimizing using indices, error: Part::pkspec1: The expression cannot be used as a part specification

I want to use Minimize where the variables to minimize are indices pointing into an array. Here a MWE that hopefully shows what my problem is. vars = u@# & /@ Range[3]; cons = Flatten@ { Table[(u[j] != #) & /@ vars[[j + 1 ;; -1]], {j, 1, 3 - 1}], 1 vec1 = {1, 2, 3}; vec2 = {1, 2, 3}; Minimize[{Total@((vec1[[#]] - vec2[[u[#]]])^2 & /@ Range[1, 3]), cons}, vars, Integers] The error I get: Part::pkspec1: The expression u[1] cannot be used as a part specification. >> Answer Ok, it seems that one can get around Mathematica trying to evaluate vec2[[u[1]]] too early by using the function Indexed[vec2,u[1]] . The working MWE would then look like the following: vars = u@# & /@ Range[3]; cons = Flatten@{ Table[(u[j] != #) & /@ vars[[j + 1 ;; -1]], {j, 1, 3 - 1}], 1 vec1 = {1, 2, 3}; vec2 = {1, 2, 3}; NMinimize[ {Total@((vec1[[#]] - Indexed[vec2, u[#]])^2 & /@ R...
Post a Comment
Read more

plotting - Plot 4D data with color as 4th dimension

I have a list of 4D data (x position, y position, amplitude, wavelength). I want to plot x, y, and amplitude on a 3D plot and have the color of the points correspond to the wavelength. I have seen many examples using functions to define color but my wavelength cannot be expressed by an analytic function. Is there a simple way to do this? Answer Here a another possible way to visualize 4D data: data = Flatten[Table[{x, y, x^2 + y^2, Sin[x - y]}, {x, -Pi, Pi,Pi/10}, {y,-Pi,Pi, Pi/10}], 1]; You can use the function Point along with VertexColors . Now the points are places using the first three elements and the color is determined by the fourth. In this case I used Hue, but you can use whatever you prefer. Graphics3D[ Point[data[[All, 1 ;; 3]], VertexColors -> Hue /@ data[[All, 4]]], Axes -> True, BoxRatios -> {1, 1, 1/GoldenRatio}]
Post a Comment
Read more