tf.math.bessel_i0 overflows to inf at 713.0 although the float64 result is finite

Issue type

Bug

Have you reproduced the bug with TensorFlow Nightly?

Yes

Source

source

TensorFlow version

tf 2.21.0

Custom code

Yes

OS platform and distribution

Linux Ubuntu 22.04

Mobile device

No response

Python version

Python 3.11

Bazel version

No response

GCC/compiler version

No response

CUDA/cuDNN version

N/A

GPU model and memory

N/A

Current behavior?

tf.math.bessel_i0 overflows to inf for a float64 input whose correct result is still finite.

For input 713.0, TensorFlow returns inf. A high-precision mpmath reference rounded to float64 is approximately 6.705128263670996e+307, which is finite and within the float64 range. This appears to be a premature overflow in tf.math.bessel_i0.

Expected behavior?

tf.math.bessel_i0(713.0) with float64 input should return a finite value close to the high-precision reference, approximately 6.705128263670996e+307, instead of overflowing to inf.

Standalone code to reproduce the issue

import os
os.environ.setdefault("TF_CPP_MIN_LOG_LEVEL", "3")
os.environ["CUDA_VISIBLE_DEVICES"] = "-1"

import numpy as np
import tensorflow as tf

x = tf.constant([713.0], dtype=tf.float64)
y = tf.math.bessel_i0(x).numpy()

truth = np.array([6.705128263670996e307], dtype=np.float64)

print("input:", x.numpy())
print("tensorflow:", y)
print("truth:", truth)
print("isfinite tensorflow:", np.isfinite(y))
print("isfinite truth:", np.isfinite(truth))

if np.isinf(y[0]) and np.isfinite(truth[0]):
    print("BUG REPRODUCED: tf.math.bessel_i0 overflows although the float64 result is finite")
else:
    print("not reproduced")


Relevant log output


input: [713.]
tensorflow: [inf]
truth: [6.70512826e+307]
isfinite tensorflow: [False]
isfinite truth: [ True]
BUG REPRODUCED: tf.math.bessel_i0 overflows although the float64 result is finite

Relevant log output

input: [713.]
tensorflow: [inf]
truth: [6.70512826e+307]
isfinite tensorflow: [False]
isfinite truth: [ True]
BUG REPRODUCED: tf.math.bessel_i0 overflows although the float64 result is finite