Skip to content

Always use vectorized kernels for exp2/exp10/log2/log10/pow/sin/cos/tan - #109717

Open
raimannma wants to merge 24 commits into
ClickHouse:masterfrom
raimannma:simd-float-math
Open

Always use vectorized kernels for exp2/exp10/log2/log10/pow/sin/cos/tan#109717
raimannma wants to merge 24 commits into
ClickHouse:masterfrom
raimannma:simd-float-math

Conversation

@raimannma

@raimannma raimannma commented Jul 7, 2026

Copy link
Copy Markdown
Contributor

Motivation

We evaluate a whole class of Float64 math functions with scalar libm, one element at a time, even though a sibling set of functions in the same binary already runs through the vectorized FastOps path. Benchmarking each in isolation (SELECT sum(fn(x)) FROM numbers(100e6), single core) two sharp tiers appear:

  • Already vectorized (FastOps): sqrt, sigmoid, tanh, exp, log — fast.
  • Still scalar libm: log2, exp2, log10, exp10, pow (and the transcendentals) — 3-12x slower.

The split isn't explained by mathematical difficulty — exp2 is inherently cheaper than exp yet runs several times slower here, and pow is the slowest of all. The scalar tier is leaving a large, consistent amount of CPU on the floor.

Conveniently, exp2/exp10 already exist in our FastOps fork (contrib/fastops) but were never wired up on the ClickHouse side; log2/log10 are just ln(x) rescaled by a constant; and pow(const_base, x) is exp2(x·log2(base)). No contrib/fastops changes are needed.

What this PR does

Routes these functions through vectorized kernels unconditionally, the same way exp/log already use FastOps with no setting (the precise scalar libm path remains only for builds without USE_FASTOPS):

  • exp2, exp10NFastOps::Exp2 / NFastOps::Exp10 (already in the fork);
  • log2, log10NFastOps::Log rescaled by 1/ln(2) / 1/ln(10);
  • pow → a fast path for the two common special cases:
    • constant positive base b: b^y = exp2(y·log2(b)) (the ~12x case);
    • constant integer exponent n with |n| <= 64: exponentiation by squaring (pow(x, 2) becomes a single multiply);
    • everything else (e.g. pow(x, 0.5), non-positive constant base, two non-constant columns, and pow(const, const) which is folded through the default constant implementation) uses precise std::pow and is bit-identical to before;
  • sin, cos, tan → a branch-free, auto-vectorized polynomial kernel (src/Functions/FastTrig.h): three-term Cody-Waite reduction to a quadrant of pi/2 plus the Cephes minimax polynomials on [-pi/4, pi/4]. No intrinsics and no runtime dispatch, so clang vectorizes it on every target. Arguments with |x| > 1e8, NaN and Inf are recomputed with libm. This addresses sin and cos call libm once per row — a vectorized kernel behind an opt-in setting is 1.4x–2.1x faster #116082 (contrib/fastops has no Sin/Cos, hence the hand-written kernel).

The unary functions plug into the existing FunctionMathUnary template via a small VectorizedFloat64Impl<Name, Kernel> adapter, so the shared template used by ~20 functions is untouched.

Benchmarks

SELECT sum(fn(x)) FROM numbers(100e6), single core, x86-64-v3 release build, AMD Ryzen 9 3900X (min-of-3, seconds):

function precise (libm) vectorized speedup
exp2 0.460 0.139 3.3x
exp10 0.985 0.144 6.8x
log2 0.522 0.208 2.5x
log10 0.464 0.192 2.4x
sin 1.074 0.329 3.3x
cos 1.076 0.323 3.3x
tan 1.094 0.344 3.2x
pow(2, x) (const base) 1.879 0.158 11.9x
pow(x, 2) (const int exponent) 0.514 0.086 6.0x
pow(x, 17) (const int exponent) 1.944 0.400 4.9x
pow(x, 0.5) (fallback) 0.604 0.587 1.0x (unchanged)
pow(x, y) (two columns, fallback) 1.854 1.851 1.0x (unchanged)
sin, |x| > 1e8 (libm fallback) 1.452 1.739 0.84x
exp (already vectorized, reference) 0.141 0.145

The sin x10 nested query from tests/performance/vectorized_trig.xml (3M rows, max_threads = 1): 0.279 s → 0.086 s (3.2x); the issue reports 1.4-2.1x on aarch64. The only slowdown is the |x| > 1e8 trig fallback, which pays for the vectorized pass before recomputing with libm.

Accuracy

Vectorized vs. precise (libm), measured in-engine over the finite domain (10M points per function, wide ranges):

function median rel. err worst rel. err notes
exp2 6.3e-13 8.9e-13 bit-exact on integer inputs (exp2(10) = 1024)
exp10 6.2e-13 9.6e-13 exp10(2) is now 99.9999999999769
log2 1.2e-12 4.6e-9 matches the accuracy of the already-shipping log/ln (it is that kernel)
log10 1.2e-12 4.6e-9
pow(b, x) const base 6.4e-13 9.1e-13
pow(x, 2) 0 0 (bit-exact)
pow(x, n) const int exponent, |n| <= 64 3e-16 - 2.5e-15 1.3e-14 repeated multiplication; the error grows like |n| * eps, which is what bounds the path at 64
sin, cos (|x| <= 1e8) 0 4.3e-16 (≤ 2 ulp) absolute error <= 1.2e-16 over the whole range; 83% of results bit-exact
tan (|x| <= 1e8) 0 5.7e-16 (≤ 2 ulp)
sin, cos, tan (|x| > 1e8, NaN, Inf) 0 (bit-exact) precise libm fallback
pow(x, y) other 0 (bit-exact) precise fallback: non-integer constant exponent, |n| > 64, non-positive constant base, two non-constant columns, or both constant

All non-finite results (overflow/underflow boundaries, NaN, ±Inf, -0) agree exactly with libm in every case; only the last few mantissa bits of finite results differ. Tests that used exp10 only as a power-of-ten generator (03022, 03268) now parse 1eN literals, 00534_exp10 and 00536_int_exp check against a 1e-11 tolerance, and a new test 04508_vectorized_float_math covers the exact cases, special values, Float32/integer inputs, and the fallbacks.

Closes: #116082

Changelog category (leave one):

  • Performance Improvement

Changelog entry (a user-readable short description of the changes that goes into CHANGELOG.md):

exp2, exp10, log2, log10 and pow are now evaluated through the vectorized FastOps kernels instead of scalar libm (2.5-12x faster; relative error ~1e-12 for exp2/exp10/pow, ~5e-9 for log2/log10, the same as the already-vectorized exp/log), and sin, cos, tan through an auto-vectorized polynomial kernel that is ~3x faster and accurate to 1-2 ulp for |x| <= 1e8. Results of these functions may differ from previous versions in the last mantissa bits (e.g. exp10(2) returns 99.9999999999769).

Documentation entry for user-facing changes

  • Documentation is written

Workflow [PR]
Sync PR [sync-upstream/pr/109717]

Route exp2, exp10, log2, log10 and pow through the vectorized FastOps
kernels (the same ones already used by exp/log) when the new
`fast_float_math` setting is enabled: 2.5-11x faster at a relative error
of ~1e-12 (exp2/exp10) to ~5e-9 (log2/log10). For pow, a constant positive
base uses exp2(y*log2(b)) and a constant integer exponent uses
exponentiation by squaring; other cases fall back to precise pow.

Default false keeps the precise scalar libm results bit-for-bit, so
existing tests are unchanged. No contrib/fastops changes are required.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
@murphy-4o murphy-4o added the can be tested Allows running workflows for external contributors label Jul 13, 2026
@clickhouse-gh

clickhouse-gh Bot commented Jul 13, 2026

Copy link
Copy Markdown
Contributor

Workflow [PR], commit [74c0308]

Summary:

job_name test_name status info comment
Fast test FAIL
04092_kql_conformance_baby_kusto FAIL cidb
00536_int_exp FAIL cidb
04093_kql_conformance_kusto_loco FAIL cidb
Code Review DROPPED
Fast test (arm_darwin) DROPPED
Build (amd_debug) DROPPED
Build (amd_asan_ubsan) DROPPED
Build (amd_tsan) DROPPED
Build (amd_msan) DROPPED
Build (amd_binary) DROPPED
Build (arm_debug) DROPPED
Build (arm_asan_ubsan) DROPPED

@clickhouse-gh clickhouse-gh Bot added the pr-feature Pull request with new product feature label Jul 13, 2026
raimannma and others added 2 commits July 13, 2026 12:10
The arm_tidy CI build failed with modernize-use-std-numbers on the
hand-written 1/ln(2) and 1/ln(10) literals. Replace them with
std::numbers::log2e and std::numbers::log10e, which are identical.

https://s3.amazonaws.com/clickhouse-test-reports/json.html?PR=109717&sha=3e84532648a56fa986924698c05170f5d4f6540f&name_0=PR&name_1=Build%20%28arm_tidy%29
ClickHouse#109717

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Comment thread src/Functions/pow.cpp
Comment thread src/Functions/pow.cpp Outdated
The `fast_float_math` pow fast path had two issues found in review:

- The negative integer exponent case computed `x^|n|` first and inverted
  at the end, so the intermediate could overflow to `+Inf` and `1/Inf`
  collapsed representable subnormals to `0` (e.g. `pow(65698.5552524023369, -64)`
  returned `0` instead of `4.74709109243818793e-309`). Reciprocate the base
  up front instead (`x^-n = (1/x)^n`) so the powers underflow gracefully.

- `FunctionPowFast` no longer inherits `FunctionMathBinaryFloat64` and lost
  its `getReturnTypeForDefaultImplementationForDynamic` override, so
  `pow(Dynamic, ...)` returned `Dynamic` instead of `Nullable(Float64)`.
  Restore the override.

Add regression tests for both around the underflow boundary and the
Dynamic return type.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Comment thread src/Functions/pow.cpp
raimannma and others added 2 commits July 13, 2026 14:37
The `fast_float_math` integer-exponent fast path in `pow` computes the power by
repeated multiplication (exponentiation by squaring), which is accurate but not
bit-identical to precise `std::pow`. Over finite integer exponents up to
magnitude 64, ~89% of results differ from `std::pow`, by up to ~100 ULP
(~1e-14 relative). For example, `pow(-0.8157093076673938, 17)` returns
`-0.031340224536743344` on the fast path versus `-0.031340224536743337` from
precise `pow`.

The setting documentation claimed `pow` is "exact for integer exponents", which
overstated the guarantee and contradicted the following sentence in the same
block ("only the last few mantissa bits of finite results differ"). Since
`fast_float_math` explicitly opts into approximate results, the fix is to state
the real accuracy rather than route integer exponents back through libm (which
would defeat the point of the fast multiplication path).

- Settings.cpp: state ~`1e-14` for integer exponents instead of "exact".
- 04508_fast_float_math: pin the real contract - agreement with precise `pow`
  to ~1e-13 relative over nontrivial floating bases, including the reported
  case. The precise reference is obtained by making the exponent a non-constant
  column (`materialize`), routing `pow` through the `std::pow` fallback even
  under `fast_float_math = 1`.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
The test asserted that the `fast_float_math` `pow` result for the specific
reviewed value `pow(-0.8157093076673938, 17)` is not bit-identical to precise
`std::pow`. Whether exponentiation-by-squaring coincides bit-for-bit with the
platform's `libm` `pow` is architecture/`libm`-dependent: on the CI runner the
two came out equal, so the `!=` half of the conjunction returned `0` and the
test failed.

Drop the fragile inequality check and keep only the real contract: agreement
with precise `pow` to ~1e-13 relative.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Comment thread src/Functions/pow.cpp Outdated
Comment thread src/Core/Settings.cpp Outdated
Comment thread src/Core/SettingsChangesHistory.cpp Outdated
The `SettingsChangesHistory.cpp` entry claimed `pow` is exact for integer
exponents, but the fast path computes them by repeated multiplication, which is
accurate to ~1e-14 and not bit-identical to precise `pow`. That is what
`Settings.cpp` documents and what `04508_fast_float_math.sql` pins.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
Comment thread src/Functions/pow.cpp Outdated
The bound was an undocumented magic number: `fast_float_math` promised the
repeated-multiplication path for "integer exponents", but `pow(x, 65)` silently
fell back to precise `pow`.

Keep the bound - it is load-bearing on two counts:

- Accuracy. Exponentiation by squaring squares its running base, so each
  squaring doubles the error already accumulated in it and the relative error
  decays like `|n| * eps`, not `log2(n) * eps`. Measured worst case over bases
  in [0.5, 2]: 6.6e-15 at |n| = 64, 1.3e-14 at 128, 1.0e-13 at 1024. 64 is where
  the ~1e-14 documented for the setting stops holding.
- Well-definedness. `y != floor(y)` does not reject integral-but-huge values:
  1e19 passes it and overflows the `Int64` cast. The bound is the only guard.

Name the constant, record the rationale, spell the limit out in the setting
text, and pin both properties with tests: `pow(x, 65)` against a non-constant
base must stay bit-identical to precise `pow` (this fails if the bound is
widened), and an exponent too large for `Int64` must not reach the integer path.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
@raimannma

raimannma commented Jul 14, 2026

Copy link
Copy Markdown
Contributor Author

Documented and tested it instead of removing it — the bound is load-bearing:

  • Accuracy. Squaring doubles the running base's accumulated error, so it decays like |n| * eps, not log2(n) * eps. Measured vs std::pow: 1.7e-15 at |n| = 17, 6.6e-15 at 64, 1.3e-14 at 128, 1.0e-13 at 1024. The setting documents ~1e-14, so 64 is where that claim stops holding.
  • Int64 safety. y != floor(y) doesn't reject 1e19 — it's integral, passes, and overflows the cast. The bound is the only guard.

So: named the constant with the rationale, spelled the limit out in the fast_float_math text, and pinned both in 04508_fast_float_math. Verified the pow(x, 65) test isn't tautological by widening the bound to 128 and confirming it fails.

Also caught: the PR description claimed integer exponents were bit-exact vs libm — only true for n = 0, 1, 2. Fixed.

64 is still a judgment call — 1024 would sit at ~1e-13, better than the ~1e-12 the const-base path already admits. Happy to raise it.

Comment thread src/Core/SettingsChangesHistory.cpp Outdated
raimannma and others added 3 commits July 15, 2026 12:04
The history entry promised ~1e-14 for integer exponents in general, but the
repeated-multiplication path only covers a constant integer exponent with
|n| <= 64; larger exponents fall back to precise pow. Settings.cpp and
04508_fast_float_math.sql already reflect that.

Co-Authored-By: Claude Opus 4.8 <noreply@anthropic.com>
raimannma and others added 3 commits July 26, 2026 17:28
Master bumped the version to 26.8, so the entry for the new
`fast_float_math` setting belongs in the 26.8 section instead of 26.7.

Co-Authored-By: Claude Opus 5 <noreply@anthropic.com>
@clickhouse-gh

clickhouse-gh Bot commented Jul 30, 2026

Copy link
Copy Markdown
Contributor

LLVM Coverage Report

Metric Baseline Current Δ
Lines 86.40% 86.50% +0.10%
Functions 91.80% 91.90% +0.10%
Branches 78.50% 78.70% +0.20%

Changed lines: Changed C/C++ lines covered: 121/123 (98.37%) · Uncovered code

Full report · Diff report

Comment thread src/Core/Settings.cpp Outdated
DECLARE(Bool, decimal_check_overflow, true, R"(
Check overflow of decimal arithmetic/comparison operations
)", 0) \
DECLARE(Bool, fast_float_math, false, R"(

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

fast_float_math changes query semantics on remote shards, but it is still declared with flags 0 instead of IMPORTANT. BaseSettings::write/read only turn an unknown setting into an exception when that bit is set; otherwise an older server in a mixed-version Distributed query just warns and ignores it.

That means one query can silently mix precise exp2/exp10/log2/log10/pow results from old shards with approximate results from new ones. Please mark this setting IMPORTANT so older nodes reject the query instead of running with divergent semantics.

Comment thread src/Functions/pow.cpp
ColumnPtr base_col = castColumn(arguments[0], f64);
ColumnPtr exp_col = castColumn(arguments[1], f64);

if (isColumnConst(*exp_col))

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This makes pow depend on ColumnConst status rather than only on argument values. With fast_float_math = 1, pow(b, 17) takes the repeated-multiplication path here, while pow(b, materialize(17)) falls back to precise std::pow; 04508_fast_float_math.sql already has to allow a tolerance between the two.

That breaks the IFunction::isDeterministic contract: constness is not stable across plan and pipeline boundaries, so a harmless rewrite/materialization can now change query results. Can we keep one math kernel per semantic case here, instead of choosing different results based on whether the exponent/base happens to be stored as ColumnConst?

raimannma and others added 4 commits August 28, 2026 10:44
`sin`, `cos` and `tan` were the remaining math functions in `expression_actions` that
evaluate one out-of-line libm call per row. With `fast_float_math = 1` they now go through
an auto-vectorized, branch-free kernel (`FastTrig.h`): a three-term Cody-Waite reduction to
a quadrant of pi/2 plus the Cephes minimax polynomials on [-pi/4, pi/4]. No intrinsics and
no runtime dispatch; clang vectorizes the loop on every target. Arguments with |x| > 1e8,
NaN and Inf are recomputed with libm, so only the finite fast range is approximate:
~1 ulp for `sin`/`cos`, ~2 ulp for `tan`, with the absolute error below 1.2e-16 everywhere.

`FastMathUnaryImpl` and `createGatedMathUnary` no longer depend on `USE_FASTOPS`, as the
trigonometric kernel does not use FastOps.

Closes: ClickHouse#116082

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01LijveiNzayj22LjT1DQGQD
@raimannma raimannma changed the title Add fast_float_math setting for vectorized exp2/exp10/log2/log10/pow Add fast_float_math setting for vectorized exp2/exp10/log2/log10/pow/sin/cos/tan Aug 28, 2026
@clickhouse-gh

clickhouse-gh Bot commented Aug 28, 2026

Copy link
Copy Markdown
Contributor

Build profile diff (arm_release)

No arm_release build profile data for commit 74c0308 - the build was skipped, reused from cache, or predates profile upload.

@CLAassistant

CLAassistant commented Aug 29, 2026

Copy link
Copy Markdown

CLA assistant check
All committers have signed the CLA.

@raimannma raimannma changed the title Add fast_float_math setting for vectorized exp2/exp10/log2/log10/pow/sin/cos/tan Always use vectorized kernels for exp2/exp10/log2/log10/pow/sin/cos/tan Aug 29, 2026
@clickhouse-gh clickhouse-gh Bot added pr-performance Pull request with some performance improvements and removed pr-feature Pull request with new product feature labels Aug 29, 2026
`exp2`, `exp10`, `log2`, `log10`, `pow`, `sin`, `cos` and `tan` now always use
the vectorized kernels (`FastOps` / `FastTrig`), the same way `exp` and `log`
already do; the precise scalar `libm` path remains only for builds without
`USE_FASTOPS`. This drops the create-time dispatcher (`createGatedMathUnary`,
`fastFloatMathEnabled`, `FunctionMathUnary.cpp`) in favor of a single
`VectorizedFloat64Impl<Name, Kernel>` template.

`pow` now uses the default implementation for constants, so
`pow(const, const)` folds to a constant column again (needed by
`numbers(pow(2, 32) - 8, 8)` etc.) and is evaluated with precise `std::pow`.

Tests that used `exp10` only as a power-of-ten generator now parse `1eN`
literals instead; `00534_exp10` and `00536_int_exp` check against a
tolerance; references with last-digit differences are updated.

Co-Authored-By: Claude Fable 5 <noreply@anthropic.com>
Claude-Session: https://claude.ai/code/session_01NuZ1m9mJBbzyxiGNWycByN
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

can be tested Allows running workflows for external contributors pr-performance Pull request with some performance improvements

Projects

None yet

Development

Successfully merging this pull request may close these issues.

sin and cos call libm once per row — a vectorized kernel behind an opt-in setting is 1.4x–2.1x faster

4 participants