5
$\begingroup$

Define a function F as follows:

F[α_, x_] = 
  N[Exp[x/2] Abs[Gamma[α/2 + I x /(2 Pi)]]^2, 10];

It is a product of two terms. As $x\to\infty$, the first term exponentially increases while the second term exponentially decreases. Their product moderately increases by a power-law.

If I evaluate F[6, 10000], the answer is 6.4163*10^16, which seems reasonable. However, if I evaluate F[6, 10000.1], the answer is 0. This is not expected since F is a continuous function.

I think the issue is on the precision, so I modified the definition F by removing N.

F[α_, x_] = 
  Exp[x/2] Abs[Gamma[α/2 + I x /(2 Pi)]]^2;

Then, even Gfintemp[6, 10000] // N gives the wrong answer 0..

How can I fix the problem, so that F gives a reasonable output for moderate values of $x$, written in exact number or numeric?

$\endgroup$
5
  • 2
    $\begingroup$ Try higher precision input than machine precision: F[6, 10000.1`20] $\endgroup$
    – Michael E2
    Jan 23, 2023 at 11:58
  • 1
    $\begingroup$ Set the precision inside F[] if you want to handle machine-precision inpuit. For instance, F[\[Alpha]_, x_] := With[{\[Alpha]0 = SetPrecision[\[Alpha], Infinity], x0 = SetPrecision[x, Infinity]}, N[Exp[x0/2] Abs[Gamma[\[Alpha]0/2 + I x0/(2 Pi)]]^2, 10] ]; -- Insteaed of SetPrecision could also use Rationalize which tends to change the value of the input slightly but also produce smaller denominators. $\endgroup$
    – Michael E2
    Jan 23, 2023 at 12:04
  • $\begingroup$ @MichaelE2 Thanks for your comment. However, your suggested code on F does not work; F[6, 10000.1] gives unwanted 0. $\endgroup$
    – Laplacian
    Jan 23, 2023 at 12:21
  • 1
    $\begingroup$ Works for me (see answer). Perhaps you did not copy the full code? Note that := in mine vs. = in F[] above is significant. $\endgroup$
    – Michael E2
    Jan 23, 2023 at 12:27
  • $\begingroup$ @MichaelE2 I made a mistake when copy-pasting your code. Now it works for me! Thanks so much for your comments and answer, as always! $\endgroup$
    – Laplacian
    Jan 23, 2023 at 12:30

1 Answer 1

6
$\begingroup$

Try higher precision input than machine precision. Set the precision inside the function if you want to handle machine-precision input. You can use SetPrecision[] or Rationalize[]. Rationalize[] tends to change the value of the input slightly but also produce smaller denominators.

f[\[Alpha]_, x_] := 
  With[{\[Alpha]0 = SetPrecision[\[Alpha], Infinity], 
    x0 = SetPrecision[x, Infinity]},
   N[Exp[x0/2] Abs[Gamma[\[Alpha]0/2 + I x0/(2 Pi)]]^2, 10]
   ];

f2[\[Alpha]_, x_] := 
  With[{\[Alpha]0 = Rationalize[Rationalize@\[Alpha], 0], 
    x0 = Rationalize[Rationalize@x, 0]},
   N[Exp[x0/2] Abs[Gamma[\[Alpha]0/2 + I x0/(2 Pi)]]^2, 10]
   ];

Evaluation is nonzero:

f[6, 10000.1]

(*  6.416572393*10^16  *)

Relative error between the methods on this expression is negligible at 10-digit precision:

(f2[6, 10000.1] - f[6, 10000.1])/f[6, 10000.1]

(*  0.*10^-10  *)

The difference in the input:

With[{x = 10000.1}, 
 N[Rationalize[Rationalize@x, 0] - SetPrecision[x, Infinity], 10]]

(*  -3.637978807*10^-13  *)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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