# Vastly incorrect answers obtained by increasing WorkingPrecision with modified Bessel functions

Bug introduced in 7.0 or earlier and fixed in 11.0

This is a follow-up to this question regarding numerical instabilities occurring with modified Bessel functions. In trying to explore J.M.'s answer to that question in greater depth by setting WorkingPrecision -> 50 (as advised in the other answer by xslittlegrass), I obtained vastly different (and incorrect!) results than when not explicitly setting a non-default WorkingPrecision (see below for plots). I am trying to figure out why this is occurring and was advised to make a separate question regarding this.

The function that I am attempting to plot (for fixed $b$, as a function of $a$) is

f[a_, b_] := (a^2 Exp[a^2/(8 b)] (BesselK[3/4, a^2/(8 b)] -
BesselK[1/4, a^2/(8 b)]))/(8 Sqrt[a b^3])


which is the analytic solution to

Integrate[x^2 Exp[-a x^2 - b x^4], {x, -∞, ∞}, Assumptions -> {a > 0, b > 0}]


I also plot the numerical solution to this,

g[a_,b_] := NIntegrate[x^2 Exp[-a x^2 - b x^4], {x, -∞, ∞}]


Without any regards to WorkingPrecision, I obtain the following plot:

amin = 10^-10; amax = 0.1; bval = 10^-3;
LogLogPlot[{f[a, bval], g[a, bval]}, {a, amin, amax}, PlotRange -> All]


Even with explicitly giving the variable/bounds higher precision,

amin = 150 10^-10; amax = 0.150; bval = 150 10^-3;


the plot remains the same as above. If, however, I add in WorkingPrecision -> 50, a plot with function values many orders of magnitude higher is obtained:

amin = 150 10^-10; amax = 0.150; bval = 150 10^-3;
LogLogPlot[{f[a, bval], g[a, bval]}, {a, amin, amax}, WorkingPrecision -> 50,
PlotRange -> All]


The fact that this is occurring is worrisome, as I would generally expect the results obtained with higher WorkingPrecision to be more trustworthy; in this case, however, they are completely wrong -- here's a representative plot:

bval = 10^-3;
Plot[x^2 Exp[-a x^2 - bval x^4] /. a -> 10^-5, {x, -15, 15}]


Why is this happening? Thank you very much!

• Do you mean Integrate instead of NIntegrate in your definition of g? May 24, 2016 at 4:19
• I believe he did mean NIntegrate[], @bb; he was comparing the solution I derived in the other thread with the results of numerical integration. May 24, 2016 at 4:27
• @J.M. I understand your point. However, Assumptions is not a legal option for NIntegrate. May 24, 2016 at 4:29
• Yes, prolly an error of copying, @bb... May 24, 2016 at 4:32
• Possible duplicate of LogPlot axes labels destroyed when working in high precision May 24, 2016 at 5:52

There appears to be a bug, not in Integrate or in BesselK, but in the vertical-axis Ticks of LogLogPlot. Consider the simple case,

LogLogPlot[Exp[x], {x, 10^-10, 1}, PlotRange -> All]


as it should be. However,

LogLogPlot[Exp[x], {x, 10^-10, 1}, WorkingPrecision -> 50, PlotRange -> All]


In fact, any value of WorkingPrecision except MachinePrecision triggers this error. Even WorkingPrecision -> \$MachinePrecision produces the error.

• Oh my. That would be quite unfortunate. I will have to see if this issue arises for various functions when I'm at my computer again tomorrow morning – this is something that would somehow need to be reported to the appropriate Wolfram team, right? May 24, 2016 at 4:40
• @PhysicsCodingEnthusiast Absolutely. Will you do so, or do you wish me to? Your call. May 24, 2016 at 4:41
• I can do so – is there a good way to do this? Should I provide a link to this question? I'll wait to do it until I'm at my computer again (currently on my phone) – in the meantime, perhaps others will chime in with their own findings. Thanks! May 24, 2016 at 4:44
• @Physics, yes, please do link to this question as evidence if and when you make a report. May 24, 2016 at 4:50
• @PhysicsCodingEnthusiast Report the bug here. Do not be surprised, if the response is slow. May 24, 2016 at 4:54