1
$\begingroup$

I am trying to evaluate an Integral which has products of Interpolating functions inside, using NIntegrate command. It gives a wrong answer. The code that I used for performing this task is given below:

    f[r_] = 1 - rh^4/r^4 - ((2*B^2)/(3*r^4))*Log[r/rh];
    q[r_] = 1 - ((2*B^2)/(3*r^4))*Log[r]; 
    h[r_] = 1 + (B^2/(3*r^4))*Log[r]; 
    gtt[r_] = (-r^2)*f[r]; 
    g11[r_] = r^2*h[r]; 
    grr[r_] = 1/(r^2*f[r]); 
    Ctt[l_] = 1/(2*Sqrt[f[rb[l]]]*rb[l]); 
    Cll[l_] = (-(1/2))*Sqrt[f[rb[l]]]*rb[l]; 
    C00[l_] = (1/(8*f[rb[l]]^(3/2)*h[rb[l]]^2*rb[l]^3))*
        (r0^4*f[r0]*h[r0]*(2*h[rb[l]]*rb[l]^2*Derivative[1][f][rb[l]]^2 + 
           2*f[rb[l]]^2*(4*h[rb[l]] + rb[l]*Derivative[1][h][rb[l]]) + 
           f[rb[l]]*rb[l]*(rb[l]*Derivative[1][f][rb[l]]*Derivative[1][h][rb[l]] + 
             2*h[rb[l]]*(Derivative[1][f][rb[l]] - rb[l]*Derivative[2][f][rb[l]]))) - 
         f[rb[l]]^2*rb[l]^4*(2*h[rb[l]]*rb[l]*Derivative[1][f][rb[l]]*
            (2*h[rb[l]] + rb[l]*Derivative[1][h][rb[l]]) + 
           f[rb[l]]*(8*h[rb[l]]^2 - rb[l]^2*Derivative[1][h][rb[l]]^2 + 2*h[rb[l]]*rb[l]*
              (4*Derivative[1][h][rb[l]] + rb[l]*Derivative[2][h][rb[l]])))); 
    D2[l_] = (r0^2*Sqrt[f[r0]*h[r0]]*(f[rb[l]]*h[rb[l]]*rb[l]^5*(h[rb[l]]*Derivative[1][f][rb[l]] - 
        f[rb[l]]*Derivative[1][h][rb[l]]) + r0^4*f[r0]*h[r0]*(h[rb[l]]*rb[l]*Derivative[1][f][rb[l]] + 
        f[rb[l]]*(4*h[rb[l]] + rb[l]*Derivative[1][h][rb[l]]))))/(2*(f[rb[l]]*h[rb[l]])^(5/2)*rb[l]^7); 
    rbP[l_] = Sqrt[(-gtt[rb[l]])*g11[rb[l]] + gtt[r0]*g11[r0]]/Sqrt[(-gtt[rb[l]])*g11[rb[l]]*grr[rb[l]]]; 
    B = 0.6; rh = 1.; r0 = 1.1; a = 2.3; 
    eom = D[Derivative[1][rb][l]^2 == rbP[l]^2, l]; 
    rs = NDSolve[{eom, Derivative[1][rb][0] == rbP[0], rb[0] == 1.10001}, rb, {l, -a, a}][[1]]
    rb[l_] = Evaluate[rb[l] /. rs]
    eq = (D[Cll[l]*Derivative[1][\[Xi]][l], l] - C00[l]*\[Xi][l])/Ctt[l]; 

and then I find the interpolated eigenfunctions from NDEigensystem by:


    {v, f} = NDEigensystem[{eq, DirichletCondition[\[Xi][l] == 0, True]}, \[Xi][l], {l, -a, a}, 2]
    e0[l_] = f[[1]]
    e1[l_] = f[[2]]

to use them in NIntegrate:


    Test = NIntegrate[D2[l]*e0[l]*e1[l]^2, {l, -a, a}, MaxRecursion -> 20]

which gives the answer as 1.27837, whereas the correct answer as given in the related paper is 3.33. Am I doing a mistake somewhere or is there any error in the use of interpolating functions inside NIntegrate?

$\endgroup$
22
  • 1
    $\begingroup$ Your code doesn't run on Mathematica v12.2! $\endgroup$ Mar 11 at 8:33
  • $\begingroup$ @UlrichNeumann I tried to copy-paste the working code cell by cell into StackExchange. But when I copy-paste it back to Mathematica from here, the code doesn't work. What can be the issue? Shall I convert the copied cell to raw input form and try? $\endgroup$
    – codebpr
    Mar 11 at 12:20
  • $\begingroup$ After copying your code to stackexchange, click {} or indent every codeline by 4 spaces $\endgroup$ Mar 11 at 12:39
  • $\begingroup$ When copying the code cells back to mathematica don't include "```" $\endgroup$ Mar 11 at 12:43
  • $\begingroup$ I tried doing that. The problem is that Mathematica is taking the entire code as a single input which may be the cause of the warning messages. I am still trying to ameliorate this issue and will edit the post soon. $\endgroup$
    – codebpr
    Mar 11 at 12:48

0

Your Answer

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

Browse other questions tagged or ask your own question.