Find root of function defined via NIntegrate

Question

I have a function defined as

rpd[r_, OptionsPattern[]] :=
  Module[
   {A = OptionValue[A]},
   NIntegrate[ 
    q/(q^2 + r) E^(-q^2/2)*
     Cos[3*q] BesselJ[1, A*q], {q, 0, \[Infinity]}]];

Options[rpd] = {A -> 1};

and plotting it for r \Elem {0,0.001} I see it crosses zero for some value of r (call it r0) which depends on the parameter a<1, like this:

I'd like to plot a graph of r0(a), so I'm trying to produce a table

rpdA = Table[{a, FindRoot[rpd[r, A -> a], {r, 0.000001}]}, {a, 0.1, 1,
     0.1}];

but I keep getting instead a lot of error messages, namely NIntegrate::ncvb and NIntegrate::slwcon. Could you help me sort this out?

Answer 1

I guess this works:

• Subdivide the interval at 10 to separate the significant oscillatory part from the superexponential decay part.
• Increase WorkingPrecision to handle the round-off error from the oscillatory part
• Use the secant method in FindRoot to prevent bad choices for r in trying to numerically approximate the gradient.
• ?NumericQ protection for rpd[].

ClearAll[rpd];
rpd[r0_?NumericQ, OptionsPattern[]] := 
  Module[{A = SetPrecision[OptionValue[A], 32], 
    r = SetPrecision[r0, 32]}, 
   NIntegrate[
     q/(q^2 + r) E^(-q^2/2)*Cos[3*q] BesselJ[1, A*q], {q, 0, 
      10, \[Infinity]}, MaxRecursion -> 20, WorkingPrecision -> 32, 
     PrecisionGoal -> 6] /; NumericQ[A]];

Options[rpd] = {A -> 1};

rpdA = Table[{a, 
    FindRoot[rpd[r, A -> a], {r, 0.000001, 0.0000001}]}, {a, 0.1, 1, 
    0.1}];

rpdA
(*
{{0.1, {r -> 0.0000846907}}, {0.2, {r -> 0.0000899893}},
 {0.3, {r -> 0.0000993154}}, {0.4, {r -> 0.000113453}},
 {0.5, {r -> 0.000133586}}, {0.6, {r -> 0.000161389}},
 {0.7, {r -> 0.000199148}}, {0.8, {r -> 0.000249928}},
 {0.9, {r -> 0.000317773}}, {1., {r -> 0.000407975}}}
*)