1
$\begingroup$

I need to calculate the integral of a function of interpolating functions. Being more specific, the function I want to integrate is the following

derUp[tau_] = 1/(drad[tau] + Sqrt[drad[tau]^2 + ftau[tau]])

where drad[tau_] = D[rad[tau], tau] and ftau[tau_] = 1 - (2 mm[tau])/rad[tau]. rad[tau] and mm[tau] are the numerical solutions of a system of differential equations, and therefore defined as interpolating functions.

The conceptual problem that I seem not being able to overcome/understand is that it seems to be impossible to calculate the following thing

func = Integrate[derUp[tau],{tau, 0, tau'}]

with tau' being a variable instead of a number. The reason behind this calculation is that I need to obtain func as a function of tau.

A solution I have tried to implement is FunctionInterpolation, that seems to works very well when the function to integrate depends on two variables, as specified in other posts that I have read on this forum, but that does not seem to be applicable to the one-variable case (at least I have not figured it out).

|improve this question
$\endgroup$
  • $\begingroup$ Try func[taup_?NumericQ] := NIntegrate[derUp[tau],{tau, 0, taup}] and call func[1.]. $\endgroup$ – march Jun 28 '16 at 21:40
  • $\begingroup$ Doesn't Integrate[derUp[tau],tau] (without explicit limits) give you what you need? $\endgroup$ – mikado Jun 28 '16 at 21:48
  • $\begingroup$ @mikado, being m[tau] and rr[tau] interpolating functions, Mathematica seems to be unable to calculate the integral using Integrate[derUp[tau],tau]. $\endgroup$ – Vale Jun 28 '16 at 21:57
  • $\begingroup$ @Vale. It appears that indefinite integration can be applied to something with Head InterpolatingFunction. Can you apply FunctionInterpolation to reduce your expressions to a single InterpolatingFunction? $\endgroup$ – mikado Jun 28 '16 at 22:43
  • 1
    $\begingroup$ You can convert your integral into an equivalent ODE problem and use NDSolve[] to construct the required solution. $\endgroup$ – J. M. is away♦ Jun 28 '16 at 23:34
1
$\begingroup$

There's not enough information in the question, but this approach fits one interpretation of the description of the situation. One can do the integral in the original NDSolve calls that generated rad and mm. (This may be what J.M. was suggesting in a comment. Or he may have meant executing a separate NDSolve call. I would recommend the approach given below.)

Example: The first two lines of the ODE system are uninteresting. The third line computes the integral of r[t], r'[t], and m[t].

{sol} = NDSolve[{
   r''[t] + r'[t] - 1/3 r[t] == 0, r[0] == 1, r'[0] == -1,
   m''[t] + m'[t] - 1/4 m[t] == 0, m[0] == 1/4, m'[0] == -1,
   f'[t] == 1/(r'[t] + Sqrt[r'[t]^2 + (1 - (2 m[t])/r[t])]), f[0] == 0
   },
  {r, m, f}, {t, 0, 10}]

ListPlot[{r, m, f} /. sol]

Mathematica graphics

The integral is returned as an InterpolatingFunction:

Mathematica graphics

|improve this answer
$\endgroup$
  • $\begingroup$ Thank you Michael, your approach is very elegant and effective. I didn't specify the system where m[tau] and rad[tau] came from to avoid having a very long post. $\endgroup$ – Vale Jun 29 '16 at 22:30
  • $\begingroup$ @Vale You're welcome. :) $\endgroup$ – Michael E2 Jun 29 '16 at 22:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.