diff --git a/compiler-rt/lib/builtins/CMakeLists.txt b/compiler-rt/lib/builtins/CMakeLists.txt --- a/compiler-rt/lib/builtins/CMakeLists.txt +++ b/compiler-rt/lib/builtins/CMakeLists.txt @@ -95,6 +95,7 @@ floatuntidf.c floatuntisf.c int_util.c + log10scaledint.c lshrdi3.c lshrti3.c moddi3.c diff --git a/compiler-rt/lib/builtins/log10scaledint.c b/compiler-rt/lib/builtins/log10scaledint.c new file mode 100644 --- /dev/null +++ b/compiler-rt/lib/builtins/log10scaledint.c @@ -0,0 +1,59 @@ +//===------- lib/log10scaledint.c - Scaled Integer log10() --------*- C -*-===// +// +// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. +// See https://llvm.org/LICENSE.txt for license information. +// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception +// +//===----------------------------------------------------------------------===// +// +// This file implements log10() for scaled integers up to a scale of 64. +// This implementation is based on Clay. S. Turner's fast binary logarithm +// algorithm. If the true result of log10() cannot be precisely stored in the +// result type, the value is rounded down towards negative infinity. +// +//===----------------------------------------------------------------------===// + +#include + +#include "int_lib.h" + +#define BASE 10 + +// Returns: log10(x), rounded down towards negative infinity in the scale +// provided in the 2nd argument +// +// Assumptions: x is represents an unsigned scaled integer +// scale is at most 64. When the scale is 0, this becomes a regular +// integer log10() function. + +int64_t __log10ScaledInt(uint64_t x, uint32_t scale) { + int64_t b = UINT64_C(1) << (scale - 1); + int64_t y = 0; + uint64_t oneval = UINT64_C(1) << scale; + uint64_t baseval = ((uint64_t)BASE) << scale; + + if (x == 0) + return INT64_MIN; // represents negative infinity + + while (x < oneval) { + x *= BASE; + y -= oneval; + } + + while (x >= baseval) { + x /= BASE; + y += oneval; + } + + for (size_t i = 0; i < scale; i++) { + x = (__uint128_t)x * x >> scale; + + if (x >= baseval) { + x /= BASE; + y += b; + } + b >>= 1; + } + + return y; +} diff --git a/compiler-rt/test/builtins/Unit/log10scaledint.c b/compiler-rt/test/builtins/Unit/log10scaledint.c new file mode 100644 --- /dev/null +++ b/compiler-rt/test/builtins/Unit/log10scaledint.c @@ -0,0 +1,85 @@ +// RUN: %clang_builtins %s %librt -o %t && %run %t + +#include +#include + +#define SCALE 31 +#define SCALED_ONE (UINT64_C(1) << SCALE) +#define __log10(x) __log10ScaledInt(x, SCALE) + +#define CHECK_LOG10(VAL, EXPECTED) \ + if (__log10(VAL) != EXPECTED) { \ + printf("error in __log10ScaledInt(%llu, %u). Expected %lld, found %lld.\n", VAL, SCALE, EXPECTED, __log10(VAL)); \ + return 1; \ + } + +#define CHECK_SUCCESS(x) \ + if (x) \ + return 1; + +int64_t __log10ScaledInt(uint64_t x, uint32_t scale); + +int check_log10() { + // Powers of 10. + // These values will not have the precision necessary to represent the true + // return values of 1og10() since reciprocal powers of 10 cannot be evenly + // represented with sums of recirpocal powers of 2 (1/10^n cannot be + // represented exactly with sums of 1/2^m). This is a binary log function so + // we can only attempt to represent powers of 10 as close as possible to sums + // of powers of 2. + // For example, we cannot get -3 exactly from log10(0.001), but we can + // precisely get log10((1 << SCALE) / 1000) and get a value close to -3 << SCALE. + // The error between this close result and the true value of -3 << SCALE will + // always be at most 1 if it is rounded. + CHECK_LOG10(SCALED_ONE / 10, -2147483652); // -2147483652 / 2^SCALE => -1.000... + CHECK_LOG10(SCALED_ONE / 100, -4294967317); + CHECK_LOG10(SCALED_ONE / 1000, -6442451226); + CHECK_LOG10(SCALED_ONE / 10000, -8589936177); + + // These should produce precise whole integer results since the scaled integer + // input is not a fraction. + CHECK_LOG10(SCALED_ONE, 0); + CHECK_LOG10(10 * SCALED_ONE, SCALED_ONE); + CHECK_LOG10(100 * SCALED_ONE, 2 * SCALED_ONE); + CHECK_LOG10(1000 * SCALED_ONE, 3 * SCALED_ONE); + CHECK_LOG10(10000 * SCALED_ONE, 4 * SCALED_ONE); + + // Sums of powers of 2 + CHECK_LOG10(SCALED_ONE / 2 + SCALED_ONE / 4, -268303894); + CHECK_LOG10(SCALED_ONE / 2 + SCALED_ONE / 4 + SCALED_ONE / 8, -124536758); + CHECK_LOG10(SCALED_ONE / 2 + SCALED_ONE / 4 + SCALED_ONE / 8 + SCALED_ONE / 16, -60191226); + CHECK_LOG10(SCALED_ONE / 8 + SCALED_ONE / 16, -1561217881); + + // Integers + CHECK_LOG10(3 * SCALED_ONE, 1024610092); + CHECK_LOG10(5 * SCALED_ONE, 1501026654); + CHECK_LOG10(7 * SCALED_ONE, 1814834221); + CHECK_LOG10(11 * SCALED_ONE, 2236373762); + + CHECK_LOG10(3 * SCALED_ONE + SCALED_ONE / 2 + SCALED_ONE / 4, 1232722760); + CHECK_LOG10(5 * SCALED_ONE + SCALED_ONE / 2 + SCALED_ONE / 4 + SCALED_ONE / 8, 1651431827); + CHECK_LOG10(7 * SCALED_ONE + SCALED_ONE / 2 + SCALED_ONE / 4 + SCALED_ONE / 8 + SCALED_ONE / 16, 1932056116); + CHECK_LOG10(11 * SCALED_ONE + SCALED_ONE / 8 + SCALED_ONE / 16, 2252137072); + + // Other numbers that can't evenly be broken into limited sums of powers of 2 + CHECK_LOG10(6746518852, 1067621222); // pi + CHECK_LOG10(5837465777, 932640298); // e + CHECK_LOG10(3037000499, 323228496); // sqrt(2) + + // Zero returns min val to represent -infinity + CHECK_LOG10(0, INT64_MIN); + + // Max + CHECK_LOG10(UINT64_MAX, 21333080777); + + // Min + CHECK_LOG10(1, -20040166791); + CHECK_LOG10(2, -19393709798); + + return 0; +} + +int main() { + CHECK_SUCCESS(check_log10()); + return 0; +}