# How to use Working Precision to get good result and plot?

I have a function f(x) which I would like to see going up for very low values of x: from an asymptotic study, I expect it to go to Infinity for low x; however, I can see an elbow and many oscillations. Are they due to numerical precision? How can I fix it? The function f(x) of interest is fon(beta).

 firstk = -8;
lastk = 0;
fon = 0.05379*beta/IntegralExact;
IntegrandON =
Simplify[Abs[
1 - (1 - I)/2*Sqrt[(beta/2)]*
Cosh[(1 - I)*z*
Sqrt[(beta/2)]]/(Sinh[(1 - I)/2*Sqrt[(beta/2)]])] //
ComplexExpand, beta > 0]^2;
IntegralExact = FullSimplify[Integrate[IntegrandON, {z, -1/2, 1/2}]];

Show[Flatten[{LogLogPlot[{fon}, {beta, 10^firstk, 10^lastk},
PlotLegends -> {"On-axis"}, PlotStyle -> Black,
PlotRange -> {{10^firstk, 10^lastk}, {10^-6, 10^(8)}},
GridLines -> {{46}, All}]}], ImageSize -> Large,
Frame -> True]


• Start by replacing 0.05379 in the definition of fon with an exact or high precision value. I'd suggest starting with 5379/100000. Do the same for any other inexact quantities in your definitions, if any. – MarcoB Aug 4 '17 at 15:26

I'd suggest doing all the heavy lifting numerically when you actually assign a value to beta, i.e. during plotting. This allows you to use NIntegrate at an appropriate working precision:

firstk = -20; (*Notice the even lower limit*)
lastk = 0;

Clear[IntegrandON]
(* Avoiding simplification and expansion significantly speeds up the calculation*)
IntegrandON[beta_] =
Abs[1 - (1 - I)/2*Sqrt[(beta/2)]*
Cosh[(1 - I)*z*
Sqrt[(beta/2)]]/(Sinh[(1 - I)/2*Sqrt[(beta/2)]])]^2;

ClearAll[fon]
fon[beta_?NumericQ] := 5379/100000*beta/NIntegrate[IntegrandON[beta], {z, -1/2, 1/2}, WorkingPrecision -> 15];

LogLogPlot[
fon[beta],
{beta, 10^firstk, 10^lastk},
PlotStyle -> Black, Frame -> True,
WorkingPrecision -> 30
]


IntegrandON =
Simplify[Abs[
1 - (1 - I)/2*Sqrt[(beta/2)]*
Cosh[(1 - I)*z*
Sqrt[(beta/2)]]/(Sinh[(1 - I)/2*Sqrt[(beta/2)]])] //
ComplexExpand, beta > 0]^2;

IntegralExact = FullSimplify[Integrate[IntegrandON, {z, -1/2, 1/2}]]


firstk = -8;
lastk = 0;
fon = 5379*10^-5*beta/IntegralExact;

LogLogPlot[fon, {beta, 10^firstk, 10^lastk},
PlotLegends -> Placed[{"On-axis"}, {.25, .7}],
PlotStyle -> Black,
PlotRange -> {{10^firstk, 10^lastk}, {1, 10^8}},
GridLines -> {10^Range[-7, -1, 2], All},
Frame -> True,
ImageSize -> Large,
WorkingPrecision -> 20,
PlotPoints -> 150]