[Date Index]
[Thread Index]
[Author Index]
Re: Compilation: Avoiding inlining
*To*: mathgroup at smc.vnet.net
*Subject*: [mg121459] Re: Compilation: Avoiding inlining
*From*: Oliver Ruebenkoenig <ruebenko at wolfram.com>
*Date*: Fri, 16 Sep 2011 05:46:21 -0400 (EDT)
*Delivered-to*: l-mathgroup@mail-archive0.wolfram.com
*References*: <201109150841.EAA29262@smc.vnet.net>
On Thu, 15 Sep 2011, DmitryG wrote:
> This message refers to the extensive exchange here:
>
> http://groups.google.com/group/comp.soft-sys.math.mathematica/browse_thread/thread/4ea9938ad1d523c6/735f864eab188778?lnk=gst&q=dmitryg#735f864eab188778
>
> The problem is to solve a system of N ordinary differential equations
> (ODE) with Runge-Kutta-4 method with compilation. Here is an example
> of working program
>
> ***************************************************************************
> (* Runge-Kutta-4 routine *)
> makeCompRK[fIn_] :=
> With[{f = fIn},
> Compile[{{x0, _Real,
> 1}, {t0, _Real}, {tMax, _Real}, {n, _Integer}},
> Module[{h, K1, K2, K3, K4, SolList, x = x0, t},
> h = (tMax - t0)/n;
> SolList = Table[x0, {n + 1}];
> Do[t = t0 + k h;
> K1 = h f[t, x];
> K2 = h f[t + (1/2) h, x + (1/2) K1];
> K3 = h f[t + (1/2) h, x + (1/2) K2];
> K4 = h f[t + h, x + K3];
> x = x + (1/6) K1 + (1/3) K2 + (1/3) K3 + (1/6) K4;
> SolList[[k + 1]] = x
> , {k, 1, n}];
> SolList],(*CompilationTarget->"C"*)]];
>
> (* Definition of the equations *)
> NN = 100;
> Timing[
> RHS = Quiet[
> Table[-x[[i]] Sin[
> 0.3 t]^2/(1 + 300 Sum[x[[j]], {j, 1, NN}]^2/NN^2), {i, 1,
> NN}]];] (* Quiet to suppress complaints *)
> F = Function[{t, x},
> Evaluate[RHS]]; (* Evaluate to use actual form of RHS *)
> Print[]
>
> (* Compilation *)
> tt0 = AbsoluteTime[];
> Timing[RK4Comp = makeCompRK[F];]
> AbsoluteTime[] - tt0
> Print[];
>
> (* Setting parameters and Calculation *)
> x0 = Table[
> RandomReal[{0, 1}], {i, 1, NN}]; t0 = 0; tMax = 100; n = 500;
> tt0 = AbsoluteTime[];
> Sol = RK4Comp[x0, t0, tMax, n];
> AbsoluteTime[] - tt0
>
> Print["Compilation: ", Developer`PackedArrayQ@Sol]
>
> (* Plotting time dependences *)
> x1List = Table[{1. t0 + (tMax - t0) k/n, Sol[[k + 1, 1]]}, {k, 0,
> n}];
> x2List = Table[{1. t0 + (tMax - t0) k/n, Sol[[k + 1, 2]]}, {k, 0, n}];
> x3List = Table[{1. t0 + (tMax - t0) k/n, Sol[[k + 1, 3]]}, {k, 0, n}];
> ListLinePlot[{x1List, x2List, x3List}, PlotMarkers -> Automatic,
> PlotStyle -> {Blue, Green, Red}, PlotRange -> {0, 1}]
>
> **********************************************************************************
>
> The drawback of this solution is what Oliver Ruebenkoenig called
> "inlining". For large NN the RHS of the equations is a lot the
> explicit code generated by Mathematica. This code then gets compiled
> and the size of the compiled code strongly increases with NN (it can
> be seen by the CompilePrint[RK4Comp] command). That is, we write a lot
> of code inline instead of letting it be generated in the course of the
> execution, with the compiled code being small.
>
> Mathematica's own compiler tolerates inlining, although the large code
> eats up the memory. If we uncomment the option "CompilationTarget-
>> "C"", we will see that for large NN compilation goes for ages. That
> is, inlining has to be avoided.
>
> Daniel Lichtblau has rewritten my similar program in the vectorized
> form that avoided inlining and could be compiled to C, and was fast.
> But vectorization is not always possible or easy. Oliver wrote that
> inlining can also be avoided by defining a function and then calling
> it, but there were no examples at that time. Now I am interested in
> the cases such as the program above where vectorization is unwanted.
> How to avoid inlining by defining a function in this case? I have
> tried to replace
>
> RHS = Quiet[Table[-x[[i]] Sin[0.3 t]^2/(1 + 300 Sum[x[[j]], {j, 1,
> NN}]^2/NN^2), {i, 1, NN}]];
>
> by
>
> FRHS[i_][x_List, t_] = -x[[i]] Sin[0.3 t]^2/(1 + 100 Sum[x[[i + j]],
> {j, 1, Min[3, NN - i]}]^2);
> RHS = Quiet[Table[FRHS[i][x, t], {i, 1, NN}]];
>
> This works correctly but does not compile, and the execution time is
> longer. Something is wrong here. Oliver, Daniel, and others, please,
> help!
>
> Thank you,
>
> Dmitry
>
>
Dmitry,
here is an idea that hopefully will get you started:
ClearAll[makeCompRK]
makeCompRK[fIn_] :=
With[{f = fIn},
Compile[{{x0, _Real,
1}, {t0, _Real}, {tMax, _Real}, {n, _Integer}},
Module[{h, K1, K2, K3, K4, SolList, x = x0, t},
h = (tMax - t0)/n;
SolList = Table[x0, {n + 1}];
Do[t = t0 + k h;
K1 = h f[t, x];
K2 = h f[t + (1/2) h, x + (1/2) K1];
K3 = h f[t + (1/2) h, x + (1/2) K2];
K4 = h f[t + h, x + K3];
x = x + (1/6) K1 + (1/3) K2 + (1/3) K3 + (1/6) K4;
SolList[[k + 1]] = x, {k, 1, n}];
SolList],(*CompilationTarget->"C",*)
CompilationOptions -> {"InlineCompiledFunctions" -> True}]];
(* Have a look at InlineCompiledFunctions *)
(*Definition of the equations*)
NN = 100;
Timing[RHS =
Quiet[Table[-x[[i]] Sin[0.3 t]^2/(1 +
300 Sum[x[[j]], {j, 1, NN}]^2/NN^2), {i, 1,
NN}]];] (*Quiet to suppress complaints*)
F = Function[{t, x}, Evaluate[RHS]];
>> 0.06
(*Compilation*)
tt0 = AbsoluteTime[];
Timing[RK4Comp = makeCompRK[F];]
AbsoluteTime[] - tt0
>> 0.2
cRHS = Compile[{{t, _Real, 0}, {x, _Real,
1}}, -x*(Sin[0.3*t]^2/(1 + 300*Total[x]^2/Length[x]^2))]
(*cRHS[0.3,{1,2,3,4}//N]*)
(*Compilation*)
tt0 = AbsoluteTime[];
Timing[RK4Comp2 = makeCompRK[cRHS];]
AbsoluteTime[] - tt0
>> 0.002
(*Needs["CompiledFunctionTools`"]*)
(*CompilePrint[RK4Comp2]*)
(*
Have a look at how "InlineCompiledFunctions" -> True/False influence the
output
*)
(*Setting parameters and Calculation*)
x0 = Table[RandomReal[{0, 1}], {i, 1, NN}]; t0 = 0; tMax = 100; n = 500;
tt0 = AbsoluteTime[];
Sol1 = RK4Comp[x0, t0, tMax, n];
AbsoluteTime[] - tt0
>> 0.43
tt0 = AbsoluteTime[];
Sol2 = RK4Comp2[x0, t0, tMax, n];
AbsoluteTime[] - tt0
>> 0.007
Norm[Sol1 - Sol2]
>> 10^-14
Sol = Sol2;
Print["Compilation: ", Developer`PackedArrayQ@Sol]
Some further notes:
- I did not try the ->"C" version, but I believe both
compiles need to be compiled to C to be effective. Have a play. I'd be
interested to about the results
- I understand that you wanted an non vectorized version. That in
principal would work in the same way; you could provide an example of such
a non-vectorizable function.
- I think this is a good illustration of why compilation to C is not
always beneficial. It is only beneficial if the compile time vs the runtime
is small.
Hope this gets you started.
Oliver
Prev by Date:
**Multisets package**
Next by Date:
**Re: Compilation: Avoiding inlining**
Previous by thread:
**Compilation: Avoiding inlining**
Next by thread:
**Re: Compilation: Avoiding inlining**
| |