[libclc] [libclc] Produce subnormal results for exp/exp2 on FP32 (PR #212696)
via cfe-commits
cfe-commits at lists.llvm.org
Mon Aug 3 18:54:55 PDT 2026
https://github.com/Lurie97 updated https://github.com/llvm/llvm-project/pull/212696
>From bc1d058505a140b143c6f96eadfd7e1914f1d21d Mon Sep 17 00:00:00 2001
From: jiajia Qian <jiajia.qian at nxp.com>
Date: Thu, 30 Jul 2026 14:03:12 +0800
Subject: [PATCH] [libclc] Produce subnormal results for exp/exp2 on FP32
The FP32 __clc_exp/__clc_exp2 implementations scaled the reduced result
by 2^p using integer bit manipulation (as_int(y) + (p << 23)). That
representation cannot encode a subnormal, so inputs whose result falls
below the smallest normal were flushed to zero via "x < llim ? 0.0f",
using ln(smallest_normal) as the lower limit. On devices that report
CL_FP_DENORM this loses the correct subnormal result (e.g. erfc(), which
is built on exp(), failed CTS math_brute_force with the reference
subnormal vs a returned 0).
Scale with __clc_ldexp instead, which produces subnormals naturally, so
no separate flush/fixup path is needed. The lower saturation limit is
relaxed from ln(smallest_normal) to ln(smallest_subnormal) (2^-149); only
inputs below that genuinely underflow to zero. On flush-to-zero devices
the hardware collapses the subnormal result to zero, matching the
previous behavior.
The exponent handed to ldexp is clamped, because for extreme-magnitude
inputs the float->int reduction is out of range and an ill-formed
exponent would make ldexp return NaN. Such inputs over/underflow and are
saturated to inf/zero by the ulim/llim selects anyway, so the clamp never
affects a representable result.
Verified on i.MX95 Mali-G310 (Panfrost/Rusticl), a flush-to-zero device
(reports no FP32 denormals, FTZ on): OpenCL-CTS math_brute_force exp and
exp2 pass 6/6 sub-tests (fp32, fp32 fast-relaxed, fp16). Without the
exponent clamp, exp/exp2 regressed to NaN for large-magnitude negative
inputs on this path.
---
libclc/clc/lib/generic/math/clc_exp.cl | 2 ++
libclc/clc/lib/generic/math/clc_exp.inc | 18 +++++++++++++-----
libclc/clc/lib/generic/math/clc_exp2.cl | 2 ++
libclc/clc/lib/generic/math/clc_exp2.inc | 16 ++++++++++++----
4 files changed, 29 insertions(+), 9 deletions(-)
diff --git a/libclc/clc/lib/generic/math/clc_exp.cl b/libclc/clc/lib/generic/math/clc_exp.cl
index 1e7174c7541c5..7eb2b5599939e 100644
--- a/libclc/clc/lib/generic/math/clc_exp.cl
+++ b/libclc/clc/lib/generic/math/clc_exp.cl
@@ -11,9 +11,11 @@
#include "clc/internal/clc.h"
#include "clc/math/clc_exp_helper.h"
#include "clc/math/clc_fma.h"
+#include "clc/math/clc_ldexp.h"
#include "clc/math/clc_mad.h"
#include "clc/math/math.h"
#include "clc/relational/clc_isnan.h"
+#include "clc/shared/clc_clamp.h"
#define __CLC_BODY "clc_exp.inc"
#include "clc/math/gentype.inc"
diff --git a/libclc/clc/lib/generic/math/clc_exp.inc b/libclc/clc/lib/generic/math/clc_exp.inc
index 5057bf8034e92..249a239206c38 100644
--- a/libclc/clc/lib/generic/math/clc_exp.inc
+++ b/libclc/clc/lib/generic/math/clc_exp.inc
@@ -37,15 +37,23 @@ _CLC_OVERLOAD _CLC_DEF __CLC_GENTYPE __clc_exp(__CLC_GENTYPE x) {
__CLC_GENTYPE y = 1.0f - (((-lo) - MATH_DIVIDE(t * v, 2.0f - v)) - hi);
- // Scale by 2^p
- __CLC_GENTYPE r = __CLC_AS_GENTYPE(__CLC_AS_INTN(y) + (p << 23));
+ // Scale by 2^p. Use ldexp so results in the subnormal range are represented
+ // exactly rather than being flushed to zero (an integer as_int(y) + (p << 23)
+ // scaling cannot encode subnormals). On flush-to-zero devices the hardware
+ // collapses the subnormal result to zero, matching the previous behavior.
+ // Clamp the exponent fed to ldexp: for extreme |x| the float->int reduction
+ // above is out of range, and an ill-formed exponent would make ldexp return
+ // NaN. Any such input over/underflows and is saturated to inf/zero by the
+ // ulim/llim selects below, so the clamp never affects a representable result.
+ __CLC_GENTYPE r = __clc_ldexp(y, __clc_clamp(p, (__CLC_INTN)-151, (__CLC_INTN)129));
// ln(largest_normal) = 88.72283905206835305366
const __CLC_GENTYPE ulim = 0x1.62e430p+6f;
- // ln(smallest_normal) = -87.33654475055310898657
- const __CLC_GENTYPE llim = -0x1.5d589ep+6f;
+ // ln(smallest_subnormal) = ln(2^-149) = -103.27892990343184. Below this the
+ // result genuinely underflows to zero.
+ const __CLC_GENTYPE llim = -0x1.9d1dap+6f;
- r = x < llim ? 0.0f : r;
+ r = x < llim ? (__CLC_GENTYPE)0.0f : r;
r = x < ulim ? r : __CLC_AS_GENTYPE((__CLC_UINTN)0x7f800000);
return __clc_isnan(x) ? x : r;
}
diff --git a/libclc/clc/lib/generic/math/clc_exp2.cl b/libclc/clc/lib/generic/math/clc_exp2.cl
index 0516963ab620f..8abb966cca045 100644
--- a/libclc/clc/lib/generic/math/clc_exp2.cl
+++ b/libclc/clc/lib/generic/math/clc_exp2.cl
@@ -11,10 +11,12 @@
#include "clc/internal/clc.h"
#include "clc/math/clc_exp_helper.h"
#include "clc/math/clc_fma.h"
+#include "clc/math/clc_ldexp.h"
#include "clc/math/clc_mad.h"
#include "clc/math/clc_rint.h"
#include "clc/math/math.h"
#include "clc/relational/clc_isnan.h"
+#include "clc/shared/clc_clamp.h"
#define __CLC_BODY "clc_exp2.inc"
#include "clc/math/gentype.inc"
diff --git a/libclc/clc/lib/generic/math/clc_exp2.inc b/libclc/clc/lib/generic/math/clc_exp2.inc
index 6da361a43ed4c..8678b7510ba26 100644
--- a/libclc/clc/lib/generic/math/clc_exp2.inc
+++ b/libclc/clc/lib/generic/math/clc_exp2.inc
@@ -36,13 +36,21 @@ _CLC_OVERLOAD _CLC_DEF __CLC_GENTYPE __clc_exp2(__CLC_GENTYPE x) {
__CLC_GENTYPE y = 1.0f - (((-lo) - MATH_DIVIDE(t * v, 2.0f - v)) - hi);
- // Scale by 2^p
- __CLC_GENTYPE r = __CLC_AS_FLOATN(__CLC_AS_INTN(y) + (p << 23));
+ // Scale by 2^p. Use ldexp so results in the subnormal range are represented
+ // exactly rather than being flushed to zero (an integer as_int(y) + (p << 23)
+ // scaling cannot encode subnormals). On flush-to-zero devices the hardware
+ // collapses the subnormal result to zero, matching the previous behavior.
+ // Clamp the exponent fed to ldexp: for extreme |x| the float->int reduction
+ // above is out of range, and an ill-formed exponent would make ldexp return
+ // NaN. Any such input over/underflows and is saturated to inf/zero by the
+ // ulim/llim selects below, so the clamp never affects a representable result.
+ __CLC_GENTYPE r = __clc_ldexp(y, __clc_clamp(p, (__CLC_INTN)-151, (__CLC_INTN)129));
const __CLC_GENTYPE ulim = 128.0f;
- const __CLC_GENTYPE llim = -126.0f;
+ // Smallest subnormal is 2^-149; below this the result underflows to zero.
+ const __CLC_GENTYPE llim = -149.0f;
- r = x < llim ? 0.0f : r;
+ r = x < llim ? (__CLC_GENTYPE)0.0f : r;
r = x < ulim ? r : __CLC_AS_FLOATN((__CLC_UINTN)0x7f800000);
return __clc_isnan(x) ? x : r;
}
More information about the cfe-commits
mailing list