Group Abstract Group Abstract

Message Boards Message Boards

Appropriate method for numerical integration of InterpolatingFunction

GROUPS:
I have a number of solutions of a system of ODEs, i.e. a set of InterpolatingFunction(s). Using NIntegrate I try calculate integrals of these functions, their derivatives and products of different combination of functions and their derivatives.
Say, f1 -- InterpolatingFunction, I find f1' as f1d[t_] = D[f1, t]. The functions are quiet smooth (red -- f1, blue -- f1'):

But during numerical integration a bulk of warnings appears:
NIntegrate[f1d f1, {t,0,1}, WorkingPrecision -> 100]
NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: ............
NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 400 times. Th............

It's clear, that an appropriate method for such case shoud be used, but what exactly?
POSTED BY: Konstantin Nosov
Answer
10 months ago
One possible thing to check is whether your InterpolatingFunctions are also computed with sufficiently higher precision.
POSTED BY: Ilian Gachevski
Answer
10 months ago
Ilian, Precision[] returns infinity values (for both f and f'). When I assing infinity to WorkingPrecision, the warning appears:

NIntegrate::wprec: Value of option WorkingPrecision -> \ is not a positive machine-sized real or integer

There exists any idea?
POSTED BY: Konstantin Nosov
Answer
10 months ago
In the above example f is an InterpolatingFunction generated by NDSolve with WorkingPrecision -> 50 (f is actually v2 from the ODE system posted here).
POSTED BY: Ilian Gachevski
Answer
10 months ago
Strongly related MMa.SE thread: "How to integrate functions of linearly interpolated data?"
In short, you need an InterpolatingFunction computed with at least 50 digits of Precision and use the "InterpolationPointsSubdivision" preprocessor of NIntegrate with appropriate "MaxSubregions" suboption:
f1 = NDSolve[{y''[x] + Sin[y[x]] == 0, y[0] == 0, y[10] == 0}, y, x,
   Method -> {"Shooting",
     "StartingInitialConditions" -> {y[0] == 0, y'[0] == 2}},
   WorkingPrecision -> 50][[1, 1, 2]];

NIntegrate[f1'[t]*f1[t], {t, 0, 5}, WorkingPrecision -> 50,
Method -> {"InterpolationPointsSubdivision", "MaxSubregions" -> 1 + Length[First@f1["Coordinates"]]}]
4.7667705747559250585699457377132228611164733033539
POSTED BY: Alexey Popkov
Answer
9 months ago