Interpolating an Antiderivative

I'd like to be able to make InterpolatingFunctions for antiderivatives of functions that can't be integrated symbolically. However, the following code returns several error messages:

FunctionInterpolation[Integrate[Sqrt[1 + x^3], {x, 0, t}], {t, 0, 10}]

Here are the messages:

Thread::tdlen: "Objects of unequal length in {-1.25,-0.416667,0.416667,1.25}^{} cannot be combined. "

Thread::tdlen: "Objects of unequal length in {0.223144 +3.14159\ I,-0.875469+3.14159\ I,-0.875469,0.223144}\ {}\n cannot be combined. "

Thread::tdlen: "Objects of unequal length in {-1.25,-0.416667,0.416667,1.25}\ {} cannot be combined."

General::stop: "Further output of Thread::tlden will be suppressed during this calculation."

FunctionInterpolation::nreal: Near t = 1.25, the function did not evaluate to a real number.

FunctionInterpolation::nreal: Near t = 1.3277777777777777, the function did not evaluate to a real number.

What's going on? Is there a simple way to make this work? Changing Integrate to NIntegrate doesn't help, though the error messages are different:

NIntegrate::nlim: x = t is not a valid limit of integration.

NIntegrate::nlim: x = t is not a valid limit of integration.

NIntegrate::nlim: x = t is not a valid limit of integration.

General::stop: "Further output of NIntegrate::nlim will be suppressed during this calculation. "

• your example is indeed integrated by MMA without problems an gives: (t^3)^(1/3)*Hypergeometric2F1[-(1/2), 1/3, 4/3, -t^3]. Sep 17 '14 at 7:19
• @Dr.WolfgangHintze Yes, it has no problems with the indefinite integral. The issue seems to be how I combined it with FunctionInterpolation. Strangely, other functions with more elementary integrals, e.g. Sqrt[1+x^2], do not produce any errors. Sep 18 '14 at 1:28
• (at)Jim: no, it's the definite integral. The indefinte integral is the true antiderivative, and it looks different. But, please have a look at my solution which shows (a) that in a (natural) 2-step procedure there are no problems with calculating the interpolation of your definite integral, and (b) this definite integral is different from the antiderivative. Sep 19 '14 at 8:59
• @JimBelk, This morning when I try first example, FunctionInterpolation[Integrate[Sqrt[1 + x^3], {x, 0, t}], {t, 0, 10}] does not give the errors you indicate, but FunctionInterpolation[Evaluate@Integrate[Sqrt[1 + x^3], {x, 0, t}], {t, 0, 10}] does give them. (The reason has to do with the ConditionalExpression returned when the Integrate is evaluated without assumptions or a definite value for t.) But can you confirm the same behavior? If not, what version of Mathematica are you using? Sep 19 '14 at 10:35
• (at) Michael E2: you might wish to look at my solution for Version 8: No problems encountered in interpolating both the definite integral (Jim's) and the antiderivative. Sep 19 '14 at 19:33

Use NDSolve

antiD = NDSolveValue[{f'[x] == Sqrt[1 + x^3], f == 0}, f, {x, 0, 10}]

Example usage:

Plot[antiD[x], {x, 0, 10}] Alternatively...

This works because this function can be antidifferentiated (by Mathematica).

antiD = FunctionInterpolation[
Evaluate @ Integrate[Sqrt[1 + x^3], {x, 0, t}, Assumptions -> 0 < t < 10],
{t, 0, 10}]

or...

integral[t_?NumericQ] := NIntegrate[Sqrt[1 + x^3], {x, 0, t}];
FunctionInterpolation[integral[t], {t, 0, 10}]

FunctionInterpolation evaluated its argument on a symbolic t. The pattern test ?NumericQ prevents evaluation of NIntegrate until an actual number is substituted for t. See also What are the most common pitfalls awaiting new users?. Note also that this way does many evaluations of NIntegrate, whereas the NDSolve method does just one integration.

• That doesn't really answer the question. Why doesn't the NIntegrate form do the right thing (by substituting the appropriate numerical value for $t?$) Sep 17 '14 at 1:02
• @IgorRivin I was adding that, but the OP didn't actually ask that Q. He asked about Integrate, and how to achieve his objective, which last Q I did answer. The first Q is probably not worth going into very deeply. Evaluate gets rid of the problem. Sep 17 '14 at 1:19
• When I need to generate accurate interpolation of a function I usually use Plot or Legacy'Plot for generating the points for Interpolation. FunctionInterpolation does not make adaptive refinement. Sep 17 '14 at 4:21
• @AlexeyPopkov That's a good point and I have done that, too. But NDSolve has PrecisionGoal and AccuracyGoal. I don't believe Plot has such controls (correct me, if I'm wrong -- I know it has "Refinement" and so forth, but it's a different way to control). The above and ifn = NDSolve[{y[t] == f[x[t]], x'[t] == 1, x == 0}, y, {t, 0, 10}] are a pretty simple commands, and they tend to be the way I've done things the last few years. Sep 17 '14 at 10:37
• @AlexeyPopkov It adapts step size according to error estimates. The estimates can be fooled (a little), so it's not perfect. NDSolve sampling can be quite different from Plot. For instance NDSolveValue[{y[t] == Sin[x[t]], x'[t] == 1, x == 0}, y, {t, 0, 2 Pi}] produces an interpolation with a max error of 6*10^-8 with 162 sample points whereas interpolating the Line produced by the corresponding Plot has a max error of 4*10^-6 with 431 sample points. But I don't want to exaggerate the differences; they're rather comparable, imo. Sep 17 '14 at 13:53

There are two aspects in this problem depending on the exact definition of the task, and, as we shall see, both are completely solved by MMA without any additional facility.

The aspects are

a) calculate the interpolation of the definite integral $f=\int_0^t \sqrt{1+x^3} \, dx$

b) calculate the interpolation of the antiderivative $fad=\int\sqrt{1+x^3} \, dx$

Part a) interpolation of Jim's definite integral

There are no problems in interpolating your definite integral if you do it in two steps.

Step 1:

Calculate the definite integral (Mathematica does this exactly)

f = Integrate[Sqrt[1 + x^3], {x, 0, t}, Assumptions -> t >= -1]

(* -> (t^3)^(1/3) Hypergeometric2F1[-(1/2), 1/3, 4/3, -t^3] *)

Step 2:

Do the interpolation

fi = FunctionInterpolation[f, {t, 0, 10}]

(* -> InterpolatingFunction[{{0.,10.}},"<>"] *)

And finally check the equality of f and fi:

Plot[{f, fi[t]}, {t, 0, 10}]

(* skip the picture *)

For completeness:

$Version (* -> "8.0 for Microsoft Windows (64-bit) (October 7, 2011)" *) This procedure works fine also with many other functions, e.g. integer powers other than 3 under the sqrt. Part b) interpolation of the antiderivative The antiderivative is given by the indefinite integral : fad = Integrate[Sqrt[1 + x^3], x]$\frac{2 \left(x+x^4+(-1)^{1/6} 3^{3/4} \sqrt{-(-1)^{1/6} \left((-1)^{2/3}+x\right)} \sqrt{1+(-1)^{1/3} x+(-1)^{2/3} x^2} \text{EllipticF}\left[\text{ArcSin}\left[\frac{\sqrt{-(-1)^{5/6} (1+x)}}{3^{1/4}}\right],(-1)^{1/3}\right]\right)}{5 \sqrt{1+x^3}}\$

This quantity, besides having a small imaginary part

(* -> 2.34276 - 1.17815*10^-16 I *)

has a dicontinuity in both Re and Im:

PlotLabel -> "Re and Im of the antiderivative fad of \!$$\*SqrtBox[\(1 + \ x^3$$]\)", AxesLabel -> {"x", "Re/Im[fad]"}] Discontinuities of antiderivatives are quite common, and they have been discussed in this community several times (cf. for example Mismatch between numerical and analytic evaluation of an integral)

This could potentially cause a problem for FunctionInterpolate over an interval including the jump.

But it is not the case here !

Rather, the interpolation returns the correct complex valued function, with just announcing - as it should do - an accuracy problem connected with the jump:

FunctionInterpolation::ncvb: FunctionInterpolation failed to meet the prescribed accuracy and precision goals after 6 recursive bisections near x = {1.99219}. Continuing to refine elsewhere. >>

(* InterpolatingFunction[{{0.,5.}},"<>"] *)

Visualizing both Re and Im of the interpolation

"Re and Im of the interpolation of the \nantiderivative fadi of \\!$$\*SqrtBox[\(1 + x^3$$]\)", AxesLabel -> {"x", "Re/Im[fadi]"}] It shows agreement with the original function except, as expected, in the vicinity of the jumps.

Summarizing part b)

the interpolation of the true antiderivative is done perfectly well by Mathematica with only a justified accuracy error message due to the jump.

Conclusion

I have checked the complete code with a fresh kernel and can state that no problems are encountered in version 8 for both interpretations of the problem.

Best regards, Wolfgang