Integrating of Interpolation function

Question

I am using Interpolation to construct an InterpolatingFunction from several points. I do not need a higher order InterpolationOrder than 1.

I wonder, in this case, is the contructed InterpolatingFunction simply a piecewise linear function? My goal is to Integrate that function, but due to performance reasons, I need to get the primitive function (evaluate the non-definite integral first). But I am wondering, whether this can be done and whether the result would be piecewise quadratic function so that the operation $F(x_i)$, where $F(x)=\int f(x)\,{\rm d} x$ and $f(x)$ is the InterpolatingFunction is as fast as plugging into a piecewise defined $\tilde{F}(x) = a_ix^2+b_ix +c_i$.

So in Mathematica, the code should look like:

f=Interpolation[data, InterpolationOrder -> 1];
F[x_]:=Integral[f[x],x];
result = Exp[F[#]]&/@biglistofx;

Are there any other performance issues that I should be aware when using Interpolation and its integration?

If the Interpolation does not produce a piecewise linear function, what would be my best approach to make the above calculation as fast as possible (note that Length@biglistofx $\sim 10^7$)?

Strongly related: http://mathematica.stackexchange.com/q/18863/280 — Alexey Popkov, Mar 14 '17 at 16:45

score 6 · Accepted Answer · answered Mar 14 '17 at 16:56

6

Integrate will compute the antiderivative of an InterpolatingFunction:

data = Table[Sin[x], {x, 1., 10.}];
f = Interpolation[data, InterpolationOrder -> 1];
F[x_] = Integrate[f[x], x] (* with =, not := *)

Mathematica graphics

Plot[{f[x], F[x]}, {x, 1, 10}]

Mathematica graphics

Timing (InterpolatingFunction and Exp automatically threads over lists, so there's no need for Map; it's faster this way, too):

result = Exp[F[#]] &@ RandomReal[{1, 10}, 10^7]; // AbsoluteTiming
(*  {18.2411, Null}  *)

Length@result
(*  10000000  *)

answered Mar 14 '17 at 16:56

Michael E2

235,386
17
334
747

2

You can also use Derivative[-1][f] to generate the pure function (i.e. InterpolatingFunction without the argument x) antiderivative. – Carl Woll Mar 14 '17 at 17:06
But are the InterpolatingFunction and its integral simple functions like piece wise defined linear/quadratic functions or are they more complex Mathematica internal expressions (and thus requiring more computation time)? Or in other words, can they be displayed as a piecewise definition? – atapaka Mar 14 '17 at 18:19
1

@leosenko The interpolating function is computed via an algorithm, no doubt optimized for its task. It's about 100 times slower than a simple math function like Exp and about 6 times faster than the correspond Piecewise function. I'd be surprised (and interested) if the interpolating function can be beat. – Michael E2 Mar 14 '17 at 18:51

score 1 · Answer 2 · answered Mar 18 '17 at 17:03

Apart from Michael's and Carl's suggestions, you can also use NDSolve[] to generate the interpolant corresponding to the integral of the original InterpolatingFunction[]:

data = Table[Sin[x], {x, 1., 10.}];
f = Interpolation[data, InterpolationOrder -> 1];
F = NDSolveValue[{ff'[x] == f[x], ff[1] == 0}, ff, {x, 1, 10}];

Plot[{f[x], F[x]}, {x, 1, 10}]

plot of interpolant and its integral

Integrating of Interpolation function

2 Answers2