# MSVC is experimenting with correctly-rounded math (LLVM libc)

**URL:** <https://discourse.julialang.org/t/msvc-is-experimenting-with-correctly-rounded-math-llvm-libc/139402>\
**Category:** Offtopic\
**Tags:** math\
**Created:** [September 12, 2026, 12:12pm UTC](https://discourse.julialang.org/t/msvc-is-experimenting-with-correctly-rounded-math-llvm-libc/139402 "2026-09-12T12:12:33Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![woclass](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/woclass/32/212699_2.png) [@woclass](https://discourse.julialang.org/u/woclass)\
**Post date:** [September 12, 2026, 12:12pm UTC](https://discourse.julialang.org/t/msvc-is-experimenting-with-correctly-rounded-math-llvm-libc/139402/1 "2026-09-12T12:12:33Z")

</div>

> Since 2015, the math functions in cmath and cstdlib are supplied by the C runtime via math.h and stdlib.h, which Microsoft ships as part of the Universal C Runtime (UCRT).  
> …
> 
> In the modern day, the UCRT math functions have known mathematical inaccuracies. See [Accuracy of Mathematical Functions in Single, Double, Double Extended, and Quadruple Precision](https://members.loria.fr/PZimmermann/papers/accuracy.pdf) by Gladman, et al. for more info (it’s hot off the press). These inaccuracies stem from implementations that are quite old, dating to times when accuracy was often computationally infeasible.
> 
> We found that the [LLVM C Library](https://libc.llvm.org/) had already made incredible progress tackling the [math portions](https://libc.llvm.org/headers/math/index.html) of the C runtime. The project upholds accuracy as the primary goal, aiming for correct rounding in all rounding modes, with Gladman, et al. confirming the project’s success. By targeting mathematical accuracy, the library effectively codifies a stable interface that allows for implementation flexibility and optimization. The library is actively maintained and has a healthy community supporting it.
> 
> [MSVC C++23: constexpr cmath with LLVM Libc - C++ Team Blog](https://devblogs.microsoft.com/cppblog/msvc-c23-constexpr-cmath-with-llvm-libc/)

Also note that:

> - [LLVM-libc](https://libc.llvm.org/) provides some correctly rounded functions (see links below): [claimed accuracy of LLVM functions](https://libc.llvm.org/headers/math/index.html#higher-math-functions)
> - seven binary64 functions are integrated into the **GNU libc** up from release 2.43 (…), and the following binary32 functions are integrated into the GNU libc up from release 2.42 (…)
> - **AMD libm** (since 4.2) integrates tanhf from CORE-MATH, and uses some CORE-MATH test cases for its erf function
> - the **Intel Math Library** or IML (checked with 2026.1.0) includes some (apparently undocumented) cr\_xxx functions
> 
> – [The CORE-MATH project](https://core-math.gitlabpages.inria.fr/)

And Rust libm (compiler-builtins)

- [Port the CORE-MATH version of `cbrt`](https://github.com/rust-lang/compiler-builtins/blob/9e648612293141da499883115a8300d3258939e2/libm/CHANGELOG.md?plain=1#L54)

---

<div class="post-metadata">

**Author:** ![woclass](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/woclass/32/212699_2.png) [@woclass](https://discourse.julialang.org/u/woclass)\
**Post date:** [September 12, 2026, 12:17pm UTC](https://discourse.julialang.org/t/msvc-is-experimenting-with-correctly-rounded-math-llvm-libc/139402/2 "2026-09-12T12:17:12Z")

</div>

It seems more and more projects are moving toward correctly-rounded libm.  
Does Julia need to follow this wave?

Right now Julia’s math is ported from openlibm (and FDLIBM),  
which is not correctly rounded.

And maybe now is a good time to remove openlibm as a dependency?

---

<div class="post-metadata">

**Author:** ![Palli](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/palli/32/3380_2.png) [@Palli](https://discourse.julialang.org/u/Palli)\
**Post date:** [September 12, 2026, 3:27pm UTC](https://discourse.julialang.org/t/msvc-is-experimenting-with-correctly-rounded-math-llvm-libc/139402/3 "2026-09-12T15:27:26Z")

</div>

This is certainly something to consider, because of e.g. [https://core-math.gitlabpages.inria.fr/cbrt64.pdf](https://core-math.gitlabpages.inria.fr/cbrt64.pdf) Not even OpenLibm 0.8.7 is accurate for sin, cos (or tan, or any other trigonometry seemingly like atanh), even for Float32 (I doubt 0.8.8 from last week changes anything).

While Julia doesn’t actually use OpenLibm for much if anything anymore, I believe the Julia standard library was ported from it. I see cbrt, sin and cos are basically unchanged from 9 years ago.

Intriguing: “`pown` has **strictly mandated IEEE 754 compliance** ” (unlike powi), and was only adopted in C23, was in IEEE since 2008, optional in C. Julia needs neither, uses multiple dispatch, but likely does something closer to the less accurate powi.

Julia does have cbrt accurate, but only for Float32. LLVM libc has it for Float64 too.

cbrt calls \_approx\_cbrt and \_improve\_cbrt neither in OpenLibm. [simplify cbrt codepaths](https://github.com/JuliaLang/julia/commit/f5c46c882c6f6bca2c6aeb59d12aa21b2b56cc10)

rsqrt for Float32 is accurate (0.500 in table) in LLVM and IML 2026.1.0. The latter has rsqrt for Float64 most accurate (0.501), and atan2pi (0.528). I think Julia has no version fast or slow or rsqrt (I’ve suggested that fast approximation before).

LLVM (“LLVM libc extension”) has f16sqrt correctly rounded, something to consider for Julia.

GNU libc 2.44 has lgamma, tgamma, erf and erfc accurate also for Float64. Where do you get these and compoundn?

I believe all univarate functions are already correctly rounded in _some_ math library for Float32 (reading the table i.e. where 0.500), since can be exhaustively checked, but such is not possible for Float64. So I’m intrigued to see there the proofs of rounding.

The expections or problematic I see are atan2pi, atanpi, compoundn, erfc, lgamma, log10p1, log2p1, logp1, rootn, tgamma.

j0, and j1 and y0 and y1, are hard also apparently, IML 2026.1.0 does them best.

> **[High Performance Correctly Rounded Math Libraries for 32-bit Floating Point...](https://blog.sigplan.org/2021/08/26/high-performance-correctly-rounded-math-libraries-for-32-bit-floating-point-representations/)**
>
> Everyone uses math libraries. Surprisingly, mainstream math libraries do not produce correct results for several thousands of inputs. Developers are seldom aware of them, which affects reproducibil…

> Using this approach, we have created RLIBM-32, a library containing implementations of several elementary functions that produce correctly rounded results for all inputs for the 32-bit float type and the 32-bit [posit](https://posithub.org/docs/Posits4.pdf) type.  
> ..  
> The functions in RLIBM-32 are significantly faster than the state of the art. The above graphs show the speedup of RLIBM-32’s float functions compared to glibc’s and Intel’s libm. On average, RLIBM-32 has 1.1x and 1.2x speedup over glibc’s float and double functions. RLIBM-32 has 1.5x and 1.6x speedup over Intel’s float and double functions. Overall, RLIBM-32 not only produces correctly rounded results for all inputs but it is faster than mainstream math libraries, which have been optimized for decades.

Only LLVM libc claims pow accurate; in Float32 (and 1 ULP for FLoat64) but elsewhere I see only 0.501 elsewhere for LLVM 22.1.8.

I’m a bit confused why Julia has (the OpenLibm dependency any more and) OpenLibm\_jll. This was used for 64-bit and 32-bit, and if I recall only not used for 32-bit Windows, why I suggested dropping OpenLibm and 32-bit Windows at the time, now it has been decided to drop 32-bit Windows.

> <https://github.com/JuliaLang/julia/commit/b6e5cb5a65b37af74a3625e7652a7b30fd4820db>
>
> This adds operators \`+%\`, \`-%\`, \`\*%\`, which are equivalent to the
> non-\`%\` versio…ns, but indicate an explicit semantic expectation that
> twos completement wrapping behavior is expected and correct. As
> discussed at JuliaCon 2014 and every year since, users have often
> requested a way to opt into explicit overflow checking of arithmetic,
> whether for debugging or because they have regulatory or procedural
> requirements that expect to be able to do this. Having explicit
> operators for overflowing semantics allows use cases that depend on
> overflow behavior for correct functioning to explicitly opt-out of any
> such checking.
> 
> I want to explicitly emphasize that there are no plans to change the
> default behavior of arithmetic in Julia, neither by introducing error
> checking nor by making it undefined behavior (as in C). The general
> consensus here is that while overflow checking can be useful, and would
> be a fine default, even if hardware supported it efficiently (which it
> doesn't), the performance costs of performing the check (through
> inhibition of other optimization) is too high. In our experience it also
> tends to be relatively harmless, even if it can be a very rude awakeing
> to users coming from Python or other languages with big-default
> integers.
> 
> The idea here is simply to give users another tool in their arsenal for
> checking correctness. Think sanitizers, not language change. This PR
> includes a macro \`@Base.Experimental.make\_all\_arithmetic\_checked\`, that
> will define overrides to make arithmetic checked, but does not include
> any mechanism (e.g. #50239) to make this fast.
> 
> Co-authored-by: Keno Fischer \<Keno@users.noreply.github.com\>
> Co-authored-by: Claude Opus 4.8 \<noreply@anthropic.com\>
> Co-authored-by: OpenAI \<noreply@openai.com\>

I see Julia developers still maintain OpenLibm, for long double and more that doesn’t apply to Julia:

> <https://github.com/JuliaMath/openlibm/commit/43f291e954ae0acf6dde65fa632b392a19ab5932>
>
> powl() in ld80/e\_powl.c kept its working scratch values (z, w, W, Wa,
> Wb, ya, yb…, u) in file-scope \`static\` variables. Two threads calling
> powl() concurrently clobbered each other's scratch, so a call could
> return arbitrary garbage -- e.g. the reproducer from the issue,
> powl(0x8.779021e7c2f81b2p+14982L, -0x8.0021b03e1f821c9p-15097L), which
> is exactly 1, would instead come back as huge values or inf.
> 
> Move that scratch to local automatic variables inside powl(). The
> \`volatile\` on z is preserved (as a local volatile) because it is
> load-bearing for the rounding/underflow behaviour, exactly like the
> twom10000 constant. The math is unchanged. reducl() and powil()
> already used local automatics and were unaffected.
> 
> Add test/regression/test-222.c: it spawns 8 threads, each repeatedly
> computing a different power with a known-exact result, and fails if any
> call ever deviates. Before the fix this fails essentially every run;
> after it, it is solid. It skips (77) when long double is not the 80-bit
> type or when threads cannot be started. The regression Makefile rule
> gains -pthread (harmless for the single-threaded tests).
> 
> Co-authored-by: Claude Opus 4.8 (1M context) \<noreply@anthropic.com\>

> <https://github.com/JuliaMath/openlibm/commit/5fe399749f9276eaa0b8403e507470da05cbbb3f>

> <https://github.com/JuliaMath/openlibm/commit/9fbeafcd4f1b6ef6aa3946c1c8faead50f38a94d>
>
> The recent addition of aarch64 fenv support uses BSD-specific types
> (\_\_uint64\_t …and \_\_uint32\_t) that are not defined on Linux systems,
> causing compilation failures on aarch64-linux targets.
> 
> This commit replaces BSD-specific types with standard C99 types:
> \- \_\_uint64\_t → uint64\_t
> \- \_\_uint32\_t → uint32\_t
> 
> These standard types are available on all platforms through the
> included \<stdint.h\> header.
> 
> Additionally, this fixes the initialization of \_\_fe\_dfl\_env from
> scalar 0 to {0} to match the typedef as a non-scalar type.
> 
> Fixes build failures on aarch64-linux-gnu and other non-BSD platforms.

---

<div class="post-metadata">

**Author:** ![giordano](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/giordano/32/2166_2.png) [@giordano](https://discourse.julialang.org/u/giordano)\
**Post date:** [September 12, 2026, 3:57pm UTC](https://discourse.julialang.org/t/msvc-is-experimenting-with-correctly-rounded-math-llvm-libc/139402/4 "2026-09-12T15:57:31Z")

</div>

In case you hadn’t seen it, the paper

> **[Accuracy of Mathematical Functions in Julia](https://arxiv.org/abs/2509.05666)**
>
> Basic computer arithmetic operations, such as $+$, $\\times$, or $÷$ are correctly rounded, whilst mathematical functions such as $e^x$, $\\ln(x)$, or $\\sin(x)$ in general are not, meaning that separate implementations may provide different results...

analyses the accuracy of Julia’s mathematical functions. Apart from the hyperbolic functions which clearly need some love, almost all others have less than 1 ULP of error (the only exceptions being `exp10` in single precision with a maximum error of 1.05 ULP, and `tan` in double precision with worst error of 1.09 ULP, but note that the double precision analysis is non-exhaustive). While almost none of the functions is correctly rounded (only correctly rounded one is `sqrt`, which just calls [`llmv.sqrt`](https://llvm.org/docs/LangRef.html#llvm-sqrt-intrinsic)), I’d say they mostly do a decent job at balancing speed and accuracy.

---

<div class="post-metadata">

**Author:** ![Palli](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/palli/32/3380_2.png) [@Palli](https://discourse.julialang.org/u/Palli)\
**Post date:** [September 12, 2026, 6:23pm UTC](https://discourse.julialang.org/t/msvc-is-experimenting-with-correctly-rounded-math-llvm-libc/139402/5 "2026-09-12T18:23:24Z")

</div>

FYI: I’ve learned there is available CoreMath.jl for

> Correctly-rounded mathematical functions with CORE-MATH

Intriguingly, it’s also faster for some (Float32) values (but slower for Float64, at least corresponding), for sinh, and it and cosh and all hyperbolic are rather inaccurate in Julia:

```julia-auto
julia> a = 10.0f0;

julia> @btime cr_sinh($a)
  14.338 ns (0 allocations: 0 bytes)
11013.232f0

julia> @btime sinh($a)
  18.879 ns (0 allocations: 0 bytes)
11013.233f0

julia> a = 10.0

julia> @btime cr_sinh($a)
  29.177 ns (0 allocations: 0 bytes)
11013.232874703393

julia> @btime sinh($a)
  16.238 ns (0 allocations: 0 bytes)
11013.232874703393

```

I’m looking into `cr_cbrt` as we speak.

> julia\> time\_cbrt()
> 
> Timing cbrt over all 2^32 Float32 values…
> 
> Timing cr\_cbrt over all 2^32 Float32 values…
> 
> # Results
> 
> cbrt: 47.377 seconds  
> cr\_cbrt: 46.909 seconds  
> cbrt: 11.031 ns/call  
> cr\_cbrt: 10.922 ns/call  
> ratio: 0.990111x
> 
> Checksums:  
> cbrt: 0x00000000  
> cr\_cbrt: 0x00000000  
> (identical)  
> (0x0000000b07e441dc, 0x0000000aebf7868c)

Before other way around (might be _noisy_ machine or interference because of rand):

> Generating 1000000 random inputs…
> 
> # ====================================================================== Float32
> 
> Accuracy:  
> Julia cbrt: max error = 0 ULP  
> cr\_cbrt: max error = 0 ULP  
> Julia non-CR results: 0 / 1000000  
> cr\_cbrt non-CR results: 0 / 1000000
> 
> Speed:  
> Julia cbrt: 9.279 ns/call  
> cr\_cbrt: 9.496 ns/call  
> ratio: 1.023x

[Both numbers lower do not make sense, if I’m also testing rand with it).]

I need to check my code because before mismatch (likely only for NaNs, then that would not be worrying):

> julia\> mismatches = exhaustive\_compare()  
> Checking all 2^32 Float32 bit patterns…  
> This is 4294967296 values.
> 
> Mismatches: 8388606
> 
> First mismatch:  
> input bits: 0x7f800001  
> input: NaN  
> cbrt bits: 0x7f800001  
> cr\_cbrt bits:0x7fc00001  
> cbrt: NaN  
> cr\_cbrt: NaN  
> 0x00000000007ffffe