Skip to content
Open
Show file tree
Hide file tree
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions libc/config/linux/riscv/entrypoints.txt
Original file line number Diff line number Diff line change
Expand Up @@ -986,6 +986,7 @@ if(LIBC_COMPILER_HAS_FIXED_POINT)
libc.src.stdfix.idivulr
libc.src.stdfix.idivuk
libc.src.stdfix.idivulk
libc.src.stdfix.rdivi
)
endif()

Expand Down
1 change: 1 addition & 0 deletions libc/config/linux/x86_64/entrypoints.txt
Original file line number Diff line number Diff line change
Expand Up @@ -1019,6 +1019,7 @@ if(LIBC_COMPILER_HAS_FIXED_POINT)
libc.src.stdfix.idivulr
libc.src.stdfix.idivuk
libc.src.stdfix.idivulk
libc.src.stdfix.rdivi
)
endif()

Expand Down
8 changes: 8 additions & 0 deletions libc/include/stdfix.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -544,3 +544,11 @@ functions:
arguments:
- type: unsigned long accum
guard: LIBC_COMPILER_HAS_FIXED_POINT
- name: rdivi
standards:
- stdc_ext
return_type: fract
arguments:
- type: int
- type: int
guard: LIBC_COMPILER_HAS_FIXED_POINT
53 changes: 53 additions & 0 deletions libc/src/__support/fixed_point/fx_bits.h
Original file line number Diff line number Diff line change
Expand Up @@ -224,6 +224,59 @@ idiv(T x, T y) {
return static_cast<XType>(result);
}

LIBC_INLINE long accum nrstep(long accum d, long accum x0) {
auto v = x0 * (2.lk - (d * x0));
return v;
}

/* Divide the two integers and return a fixed_point value
*
* For reference, see:
* https://en.wikipedia.org/wiki/Division_algorithm#Newton%E2%80%93Raphson_division
* https://stackoverflow.com/a/9231996
*/
template <typename XType> LIBC_INLINE constexpr XType divi(int n, int d) {
// If the value of the second operand of the / operator is zero, the
// behavior is undefined. Ref: ISO/IEC TR 18037:2008(E) p.g. 16
LIBC_CRASH_ON_VALUE(d, 0);

if (LIBC_UNLIKELY(n == 0)) {
return FXRep<XType>::ZERO();
}
bool d_is_signed = false;
if (d < 0) {
d = (d * -1);
d_is_signed = true;
}

unsigned int nv = static_cast<unsigned int>(n);
unsigned int dv = static_cast<unsigned int>(d);
unsigned int clz = cpp::countl_zero<unsigned int>(dv) - 1;
unsigned long int scaled_val = dv << clz;
/* Scale denominator to be in the range of [0.5,1] */
FXBits<long accum> d_scaled{scaled_val};
unsigned long int scaled_val_n = nv << clz;
/* Scale the numerator as much as the denominator to maintain correctness of
* the original equation
*/
FXBits<long accum> n_scaled{scaled_val_n};
long accum n_scaled_val = n_scaled.get_val();
long accum d_scaled_val = d_scaled.get_val();
/* x0 = (48/17) - (32/17) * d_n */
long accum a = 2.8235lk; /* 48/17 */
long accum b = 1.8823lk; /* 32/17 */
long accum initial_approx = a - (b * d_scaled_val);
long accum val = nrstep(d_scaled_val, initial_approx);
val = nrstep(d_scaled_val, val);
val = nrstep(d_scaled_val, val);
val = nrstep(d_scaled_val, val);
Copy link
Contributor

Choose a reason for hiding this comment

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

Can you add error estimations for the initial approximation and after each Newton-Raphson step in the comments?

Copy link
Author

@bojle bojle Aug 30, 2025

Choose a reason for hiding this comment

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

I've added error info. I've done it to the best of my understanding, though I could use a theoretical cross-check of the explanation i've provided.

long accum res = n_scaled_val * val;
if (d_is_signed) {
res *= static_cast<XType>(-1);
}
return static_cast<XType>(res);
}

} // namespace fixed_point
} // namespace LIBC_NAMESPACE_DECL

Expand Down
14 changes: 14 additions & 0 deletions libc/src/stdfix/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -89,6 +89,20 @@ foreach(suffix IN ITEMS r lr k lk ur ulr uk ulk)
)
endforeach()

foreach(suffix IN ITEMS r)
add_entrypoint_object(
${suffix}divi
HDRS
${suffix}divi.h
SRCS
${suffix}divi.cpp
COMPILE_OPTIONS
${libc_opt_high_flag}
DEPENDS
libc.src.__support.fixed_point.fx_bits
)
endforeach()

add_entrypoint_object(
uhksqrtus
HDRS
Expand Down
21 changes: 21 additions & 0 deletions libc/src/stdfix/rdivi.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,21 @@
//===-- Implementation of rdivi function ---------------------------------===//
//
// 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
//
//===----------------------------------------------------------------------===//

#include "rdivi.h"
#include "include/llvm-libc-macros/stdfix-macros.h" // fract
#include "src/__support/common.h" // LLVM_LIBC_FUNCTION
#include "src/__support/fixed_point/fx_bits.h" // fixed_point
#include "src/__support/macros/config.h" // LIBC_NAMESPACE_DECL

namespace LIBC_NAMESPACE_DECL {

LLVM_LIBC_FUNCTION(fract, rdivi, (int a, int b)) {
return fixed_point::divi<fract>(a,b);
}

} // namespace LIBC_NAMESPACE_DECL
21 changes: 21 additions & 0 deletions libc/src/stdfix/rdivi.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,21 @@
//===-- Implementation header for rdivi ------------------------*- 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
//
//===----------------------------------------------------------------------===//

#ifndef LLVM_LIBC_SRC_STDFIX_RDIVI_H
#define LLVM_LIBC_SRC_STDFIX_RDIVI_H

#include "include/llvm-libc-macros/stdfix-macros.h"
#include "src/__support/macros/config.h"

namespace LIBC_NAMESPACE_DECL {

fract rdivi(int a, int b);

} // namespace LIBC_NAMESPACE_DECL

#endif // LLVM_LIBC_SRC_STDFIX_RDIVI_H
16 changes: 16 additions & 0 deletions libc/test/src/stdfix/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -120,6 +120,22 @@ foreach(suffix IN ITEMS r lr k lk ur ulr uk ulk)
)
endforeach()

foreach(suffix IN ITEMS r)
add_libc_test(
${suffix}divi_test
SUITE
libc-stdfix-tests
HDRS
DivITest.h
SRCS
${suffix}divi_test.cpp
DEPENDS
libc.src.stdfix.${suffix}divi
libc.src.__support.fixed_point.fx_bits
libc.hdr.signal_macros
)
endforeach()

add_libc_test(
uhksqrtus_test
SUITE
Expand Down
61 changes: 61 additions & 0 deletions libc/test/src/stdfix/DivITest.h
Original file line number Diff line number Diff line change
@@ -0,0 +1,61 @@
//===-- Utility class to test fxdivi functions ------------------*- 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
//
//===----------------------------------------------------------------------===//

#include "src/__support/CPP/type_traits.h"
#include "src/__support/fixed_point/fx_bits.h"
#include "src/__support/fixed_point/fx_rep.h"
#include "test/UnitTest/Test.h"

template <typename XType> XType get_epsilon() {
// TODO: raise error
return 0;
}
template <> fract get_epsilon() { return FRACT_EPSILON; }
template <> unsigned fract get_epsilon() { return UFRACT_EPSILON; }
template <> long fract get_epsilon() { return LFRACT_EPSILON; }

template <typename XType>
class DivITest : public LIBC_NAMESPACE::testing::Test {
using FXRep = LIBC_NAMESPACE::fixed_point::FXRep<XType>;
using FXBits = LIBC_NAMESPACE::fixed_point::FXBits<XType>;

public:
typedef XType (*DivIFunc)(int, int);

void testBasic(DivIFunc func) {
XType epsilon = get_epsilon<XType>();
EXPECT_LT((func(2, 3) - 0.666666r), epsilon);
EXPECT_LT((func(3, 4) - 0.75r), epsilon);
EXPECT_LT((func(1043, 2764) - 0.37735r), epsilon);
EXPECT_LT((func(60000, 720293) - 0.083299r), epsilon);

EXPECT_EQ(func(128, 256), 0.5r);
EXPECT_EQ(func(1, 2), 0.5r);
EXPECT_EQ(func(1, 4), 0.25r);
EXPECT_EQ(func(1, 8), 0.125r);
EXPECT_EQ(func(1, 16), 0.0625r);

EXPECT_EQ(func(-1, 2), -0.5r);
EXPECT_EQ(func(1, -4), -0.25r);
EXPECT_EQ(func(-1, 8), -0.125r);
EXPECT_EQ(func(1, -16), -0.0625r);
}

void testSpecial(DivIFunc func) {
EXPECT_EQ(func(0,10), 0.r);
EXPECT_DEATH([func] { func(10, 0); }, WITH_SIGNAL(-1));
EXPECT_EQ(func(-32768,32768), FRACT_MIN);
EXPECT_EQ(func(32767,32768), FRACT_MAX);
Comment on lines +49 to +50
Copy link
Contributor

Choose a reason for hiding this comment

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

I think this will only be true if the number of fractional bits for a fract type is 16, but it can vary depending on the platform, so let's maybe instead derive these number of fractional bits. Something like

EXPECT_EQ(func(-1 << FRACT_FBIT, 1 << FRACT_FBIT), FRACT_MIN);
EXPECT_EQ(func((1 << FRACT_FBIT) - 1, 1 << FRACT_FBIT), FRACT_MAX);

}
};

#define LIST_DIVI_TESTS(Name, XType, func) \
using LlvmLibc##Name##diviTest = DivITest<XType>; \
TEST_F(LlvmLibc##Name##diviTest, Basic) { testBasic(&func); } \
TEST_F(LlvmLibc##Name##diviTest, Special) { testSpecial(&func); } \
static_assert(true, "Require semicolon.")
14 changes: 14 additions & 0 deletions libc/test/src/stdfix/rdivi_test.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,14 @@
//===-- Unittests for rdivi -----------------------------------------------===//
//
// 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
//
//===----------------------------------------------------------------------===//

#include "DivITest.h"

#include "llvm-libc-macros/stdfix-macros.h" // fract
#include "src/stdfix/rdivi.h"

LIST_DIVI_TESTS(r, fract, LIBC_NAMESPACE::rdivi);