2
$\begingroup$

I am trying to do a numerical integration plot in Mathematica. Specifically, I want to plot $r$ versus $\tau$ where, $$\frac{dr}{d\tau}=\sqrt{2F(r)}$$ $F$ being a quartic polynomial. The limit of $r$ is supposed to be from 0 to 1, and $\tau$ is supposed to extend from something like -3 to +3. I simply cannot figure out the logic to do this. But I must perform numerical integration because $F$ is not simple enough for the integral to be solvable analytically. Please help.

$\endgroup$
2
  • $\begingroup$ Have you looked at the examples for NDSolve? $\endgroup$
    – user484
    Commented Apr 2, 2016 at 5:32
  • $\begingroup$ Actually, this DE can be solved in terms of elliptic functions, but can reduce to elementary solutions if the quartic has multiple roots. Both the Jacobi and Weierstrass functions are built-in, so this is doable in principle. The question then is what your quartic looks like in general (all real roots, one pair of complex roots, etc.) $\endgroup$ Commented Apr 2, 2016 at 8:18

1 Answer 1

1
$\begingroup$

First of all make sense of the equation, when you have a differential equation of that type you can always put it in the form $\frac{\mathrm{d}r}{\sqrt{2F\left(r\right)}}=\mathrm{d}\tau\rightarrow\int_{r_{0}}^{r}\frac{\mathrm{d}r'}{\sqrt{2F\left(r'\right)}}=\int_{\tau_0}^{\tau}\mathrm{d}\tau'$, so the idea is integrate until a certain point use the "cumulative integral" and make a list:

$list[1]=\{\int_{r_{0}}^{r_{1}}\frac{\mathrm{d}r'}{\sqrt{2F\left(r'\right)}},\int_{\tau_0}^{\tau_{1}}\mathrm{d}\tau\}$

$list[2]=\{\int_{r_{0}}^{r_{2}}\frac{\mathrm{d}r'}{\sqrt{2F\left(r'\right)}},\int_{\tau_0}^{\tau_{2}}\mathrm{d}\tau\}$,

$\hspace{3cm}\vdots$

$list[N]=\{\int_{r_{0}}^{r_{N}}\frac{\mathrm{d}r'}{\sqrt{2F\left(r'\right)}},\int_{\tau_0}^{\tau_{N}}\mathrm{d}\tau\}$,

So, for example, you want to get the integral for five values of $r$ $\rightarrow$ five values for $\tau$, thus:

F[r_] := 1/(2.0 * r^2)
tauVals = Table[i,{i,1,5}];
rVals  = Table[i,{i,0.2,1.0,0.2}]];
sols = NIntegrate[1/(Sqrt[2.0*F[r]]),{r,0,#} & /@ rVals;

Now what you want is a list plot of $\{f\left(r\right),\tau\}$, that will actally give you a quadratic funtion, in this case:

res = {{0,0},};
Do[res=Append[res,{tauVal[[i]],sols[[i]]}],{i,1,Length[tauVals]}];
ListPlot[res,AxesLabel -> {"tau","f(r)"}]

Now, if you do this with a very big number of points, you can actually interpolate the function (I only use five points as an example):

intFunc = Interpolate[res];

Now it's just a question of finding the solution to a series of points within the range and then making the plot:

rSols = r /. NSolve[intFunc[r] == #, r] & /@ tauVals;
Do[rSols[[i]] = Prepend[rSols[[i]], tauVals[[i]]], {i, 1,Length[tauVals]}];

With this interpolating functions you usually will get warings to use Reduce, but just ignore the warnings, unless you somehow run into a singularity when integrating. Plot to make sure you get a reasonable result:

Plot[Interpolation[rSols][x], {x, 1, 5}]

That is ofcourse something that looks like a square root function. I will generate the graphs later and place an update if it's the case.

$\endgroup$
1
  • $\begingroup$ Thank you. I will check and comment if it is working. $\endgroup$ Commented Apr 3, 2016 at 3:41

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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