# How to compile procedural program effectively

Consider these four function in different form:

f1[ls1_] := Append[Rest@ls1, 0.] + Prepend[Most@ls1, 0.]

f2 = Compile[{{ls1, _Real, 1}},
Append[Rest@ls1, 0.] + Prepend[Most@ls1, 0.],
CompilationTarget -> "C"];

f3[ls_] := Module[{n = Length@ls, tp = ls},
tp[[1]] = ls[[2]];
tp[[n]] = ls[[n - 1]];
Do[
tp[[i]] = ls[[i - 1]] + ls[[i + 1]], {i, 2, n - 1}
];
tp]

f4 = Compile[{{ls, _Real, 1}},
Module[{n = Length@ls, tp = ls},
tp[[1]] = ls[[2]];
tp[[n]] = ls[[n - 1]];
Do[
tp[[i]] = ls[[i - 1]] + ls[[i + 1]], {i, 2, n - 1}
];
tp], CompilationTarget -> "C"
];


Comparing their performance, we can see that the functional version f1 and f2 have about the same performance, but the compiled procedural version f4 is about 2X slower.

ls = Sin /@ (Range[1, 2000]/2000.);

AbsoluteTiming[TimeConstrained[Do[f1[ls];, {20000}], 5]]
AbsoluteTiming[TimeConstrained[Do[f2[ls];, {20000}], 5]]
AbsoluteTiming[TimeConstrained[Do[f3[ls];, {20000}], 5]]
AbsoluteTiming[TimeConstrained[Do[f4[ls];, {20000}], 5]]

(*{0.976474,Null}*)
(*{0.784508,Null}*)
(*{5.000856,$Aborted}*) (*{1.595374,Null}*)  However, if we test the procedural version in Fortran (code at the end), it runs about 4X faster than f2. xslittlegrass$ ifort -r8 -O0 test.f90
xslittlegrass$./a.out 0.232574000000000 xslittlegrass$


Questions:

1. Why f4 is about 2X slower than f1 and f2?
2. Why f4 is about 7X slower than the Fortran version?
3. Can we make a compiled version that has comparable performance to the Fortran version?
4. I thought that a procedural function compiled to C is comparable in performance compared to procedural functions written in C. But this example seems to disagree with it. So in general how large is the gap of the performance between these two?

test.f90

program main
implicit none

integer,parameter :: N0 = 2000
integer i,j
real(kind=8) :: list(N0), result(N0),start,finish

do i = 1, N0
list(i) = sin(i/real(N0))
end do

call CPU_TIME(start)

do j = 1, 20000
result(1) = list(2)
result(N0) = list(N0-1)
do i=2,N0-1
result(i) = list(i+1) + list(i-1)
end do
end do

call CPU_TIME(finish)

write(*,*) finish-start

end program

• Have you looked at the CompilePrint output for those compiled functions? Maybe something in the way they are translated to compiled code might explain the difference. Also Leonid has a good set of pointers on How to compile effectively. – MarcoB Jul 22 '15 at 21:49
• @MarcoB one looks procedural (obviously!) while the other uses compiled implementations of Inset and Drop. – dr.blochwave Jul 22 '15 at 22:02
• Are you sure that mma is only 7x slower than fortran? Which fortran compiler do you use? If you use intel cpu, you should choose intel fortran, you will be astonished by the result. On my computer, f4 takes 1.43sec, while your fortran code if compiled with gfortran, takes 0.2sec, But if compiled with ifort, only 0.0268sec !!!!! Be caution, you should add an output code to get the real time, or ifort is smart enough to know that you didn't output anything, so it will just doesn't do anything. – matheorem Jan 8 '16 at 10:42
• and I forgot to put RuntimeOption->"Speed" timing, with this option, f4 takes 0.393sec – matheorem Jan 8 '16 at 10:45
• @matheorem Thanks for the test. I was using the intel compiler, but with optimization turned off. – xslittlegrass Jan 8 '16 at 15:47

This is an answer to Q1 and partly Q4, really. I can't test your Fortran version at the moment, but it would be an interesting comparison.

You can improve the performance of f4 compared to f1 and f2 by setting RuntimeOptions -> "Speed". Clearly the change in runtime settings (mainly "CatchMachineIntegerOverflow" it seems...) from the defaults has a different effect on the two functions.

For instance:

f2 = Compile[{{ls1, _Real, 1}},
Append[Rest@ls1, 0.] + Prepend[Most@ls1, 0.],
CompilationTarget -> "C", RuntimeOptions -> "Speed"];

f4 = Compile[{{ls, _Real, 1}},
Module[{n = Length@ls, tp = ls}, tp[[1]] = ls[[2]];
tp[[n]] = ls[[n - 1]];
Do[tp[[i]] = ls[[i - 1]] + ls[[i + 1]], {i, 2, n - 1}];
tp], CompilationTarget -> "C", RuntimeOptions -> "Speed"];

AbsoluteTiming[TimeConstrained[Do[f2[ls];, {20000}], 5]]
(* 0.152 seconds *)

AbsoluteTiming[TimeConstrained[Do[f4[ls];, {20000}], 5]]
(* 0.127 seconds *)

• Hi, blochwave, if you are interested in fortran timing, see my comment to xslittlegrass – matheorem Jan 8 '16 at 10:42
• basically, if compiled with ifort on intel cpu, even RuntimeOptions->"Speed" is still more than 10x slower than fortran – matheorem Jan 8 '16 at 10:44

In the fortran version, you have the large loop within the executable and you call the executable only once. Also quite important: You allocate memory for result only once. In order to make the benchmark fair, you should compare to something like

f5 = Compile[{{ls, _Real, 1}},
Module[{n = Length@ls, tp},
tp = Table[0., {Length[ls]}];
Do[
tp[[1]] = ls[[2]];
tp[[n]] = ls[[n - 1]];
Do[tp[[i]] = ls[[i - 1]] + ls[[i + 1]], {i, 2, n - 1}];
, {20000}];
tp],
CompilationTarget -> "C",
RuntimeOptions -> "Speed"
];


On my machine, calling f5[ls] takes about 0.067 s while Do[f4[ls];, {20000}] needs about 0.58 s. The fortran variant (see below) compiled with gfortran -o bla test2.f90 -O3 needs about 0.012 s. Moreover, the call f6[ls] with the function defined below needs only 0.034 seconds which is not too far away from the fortran timing.

f6 = Compile[{{ls, _Real, 1}},
Module[{n = Length@ls, tp},
tp = Table[0., {n}];
Do[
tp[[1]] = CompileGetElement[ls, 2];
tp[[n]] = CompileGetElement[ls, n - 1];
Do[
tp[[i]] = CompileGetElement[ls, i - 1] + CompileGetElement[ls, i + 1],
{i, 2, n - 1}],
{20000}];
tp],
CompilationTarget -> "C",
RuntimeOptions -> "Speed"
];


test2.f90

program main
implicit none
integer,parameter :: N0 = 2000
integer i,j
real (kind=8) :: list(N0), result(N0),start,finish
do i = 1, N0
list(i) = sin(i/real(N0))
end do
call CPU_TIME(start)
do j = 1, 20000
result(1) = list(2)
result(N0) = list(N0-1)
do i=2,N0-1
result(i) = list(i+1) + list(i-1)
end do
end do
call CPU_TIME(finish)
write(*,*) result
write(*,*) finish-start
end program


PS.: I have still not used parallelization, here. Even with the dullest way to do that within Mathematica, I get (on a Quad Core CPU):

AbsoluteTiming[ParallelDo[f6[ls], {i, 1, $KernelCount}]][[1]]/($KernelCount)
(* 0.0125135 *)


and by even increasing the number of jobs:

AbsoluteTiming[ParallelDo[f6[ls], {i, 1, $KernelCount 10}]][[ 1]]/($KernelCount 10)
(* 0.00855335 *)


Admittedly, it is not exactly fair to compare this to the unparallelized fortran code. I added this only to show another possibility to speed up the Mathematica code.