1
$\begingroup$

I'm currently working with the HeunB function, and am finding that its numerical precision is not sufficient at large arguments - in particular I have the function $\text{HeunB}\left[-k^2,1-\beta ,\frac{1}{2},0,2 \beta ,t\right]$, which from what I can tell, for $\beta\geq 0$ and $t\gg 1/\sqrt{\beta}$, grows like $e^{-\beta t^2/2}$.

This is a problem since my calculation ultimately comes down to the cancellation between two such functions, and requires more digits than mathematica is willing to give me:

Manipulate[ LogLogPlot[{Abs[\[Xi] 2 E^(\[Beta] \[Xi]^2) ( HeunB[-k^2, 1, 3/2, 0, 2 \[Beta], t] HeunB[-k^2, 1 - \[Beta], 1/ 2, 0, 2 \[Beta], \[Xi]] - Sqrt[\[Xi]/t] HeunB[-k^2, 1, 3/2, 0, 2 \[Beta], \[Xi]] HeunB[-k^2, 1 - \[Beta], 1/2, 0, 2 \[Beta], t])]}, {\[Xi], 0, t}, PlotPoints -> 4, GridLines -> {{Sqrt[2/\[Beta]]}, {}}, PlotRange -> {0.000001, 100}], {k, 0, 10}, {\[Beta], 0.1, 10}, {t, 1, 20}]

Plotting the above function, you'll notice a dramatic numerical blow up as you take $t$ and $\beta$ larger.

Is there any way to get mathematica to spit out more digits?

$\endgroup$
3
  • 4
    $\begingroup$ With[{k = 1, β = 1/10, t = 123}, N[HeunB[-k^2, 1 - β, 1/2, 0, 2 β, t], 10^3]] correctly spits out 1000 digits of the result. What exactly is your question? Please provide a concrete example of a result you're not happy with: input, actual output, desired output. Including values of all required parameters for a concrete example. $\endgroup$
    – Roman
    May 10, 2023 at 19:11
  • $\begingroup$ See the plot I just added for a concrete example $\endgroup$
    – Guy
    May 10, 2023 at 19:16
  • 4
    $\begingroup$ You are using machine-precision numbers as arguments, and so the result will always be a machine-precision number. Try using parameters with higher precision, for example exact rationals or extended-precision numbers. $\endgroup$
    – Roman
    May 10, 2023 at 19:30

1 Answer 1

4
$\begingroup$
$Version

(* "13.2.1 for Mac OS X ARM (64-bit) (January 27, 2023)" *)

ClearAll["Global`*"];

Don't use machine precision and specify a WorkingPrecision

Manipulate[
 {k, β, t} = Rationalize[{kv, βv, tv}];
 LogLogPlot[{Abs[ξ 2 E^(β ξ^2) (HeunB[-k^2, 1, 3/2, 0, 
         2 β, t]*
         HeunB[-k^2, 1 - β, 1/2, 0, 2 β, ξ] -
       Sqrt[ξ/t] HeunB[-k^2, 1, 3/2, 0, 2 β, ξ]*
         HeunB[-k^2, 1 - β, 1/2, 0, 2 β, t])]},
  {ξ, 0, t},
  GridLines -> {{Sqrt[2/β]}, {}},
  PlotRange -> {10^-8, 100},
  MaxRecursion -> 5,
  WorkingPrecision -> 25],
 {{kv, 10, "k"}, 0, 10, 0.1, Appearance -> "Labeled"},
 {{βv, 10, "β"}, 0.1, 10, 0.1, Appearance -> "Labeled"},
 {{tv, 10, "t"}, 1, 20, 0.1, Appearance -> "Labeled"},
 SynchronousUpdating -> False,
 TrackedSymbols :> {kv, βv, tv}]

enter image description here

EDIT: To get the full plot, i.e., out to ξ = t, much higher precision is required. This significantly slows down the plotting. The plot is split in two parts so that the higher precision is only used for the larger values of ξ

Manipulate[{k, β, t} = Rationalize[{kv, βv, tv}];
 Block[{$MaxExtraPrecision = 200},
  Show[
   LogLogPlot[
    {Abs[ξ 2 E^(β ξ^2) (HeunB[-k^2, 1, 3/2, 0, 2 β, t]*
          HeunB[-k^2, 1 - β, 1/2, 0, 2 β, ξ] -
         Sqrt[ξ/t] HeunB[-k^2, 1, 3/2, 0, 2 β, ξ]*
          HeunB[-k^2, 1 - β, 1/2, 0, 2 β, t])]},
    {ξ, 0, t/4},
    GridLines -> {{Sqrt[2/β]}, {}},
    PlotRange -> {{0, t}, {10^-8, 300}},
    MaxRecursion -> 5,
    WorkingPrecision -> 25,
    Exclusions -> None],
   LogLogPlot[
    {Abs[ξ 2 E^(β ξ^2) (HeunB[-k^2, 1, 3/2, 0, 2 β, t]*
          HeunB[-k^2, 1 - β, 1/2, 0, 2 β, ξ] -
         Sqrt[ξ/t] HeunB[-k^2, 1, 3/2, 0, 2 β, ξ]*
          HeunB[-k^2, 1 - β, 1/2, 0, 2 β, t])]},
    {ξ, t/4, t},
    GridLines -> {{Sqrt[2/β]}, {}},
    PlotRange -> {{0, t}, {10^-8, 300}},
    MaxRecursion -> 5,
    WorkingPrecision -> 400,
    Exclusions -> None]]],
 {{kv, 10, "k"}, 0, 10, 0.1, Appearance -> "Labeled"},
 {{βv, 10, "β"}, 0.1, 10, 0.1, Appearance -> "Labeled"},
 {{tv, 10, "t"}, 1, 20, 0.1, Appearance -> "Labeled"},
 SynchronousUpdating -> False,
 TrackedSymbols :> {kv, βv, tv}]

enter image description here

$\endgroup$
4
  • $\begingroup$ I initially thought this was the solution, but it doesn't seem this works at larger $\xi$, (I believe your plot is cut off before $\xi = t$) $\endgroup$
    – Guy
    May 10, 2023 at 21:46
  • $\begingroup$ Oh, i suppose the resolution is just ridiculous working precision $\endgroup$
    – Guy
    May 10, 2023 at 21:52
  • $\begingroup$ No, the plot resolution is also due to elimination of PlotPoints -> 4 and addition of MaxRecursion -> 5 $\endgroup$
    – Bob Hanlon
    May 10, 2023 at 23:27
  • $\begingroup$ Oh, I meant the resolution to my problem, not of the plot lol - thanks for the help! $\endgroup$
    – Guy
    May 11, 2023 at 17:01

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