CPU backend: math functions are scalar libm calls per work-item; PoCL never uses its vector math library #860

On the PoCL-based CPU backend, sin, cos, exp, log, pow, tan, … in a kernel become one scalar libm call per work-item. They are 1.5–2.5x slower than the same loop in plain Julia, and keep PoCL from vectorizing the work-item loop. Ideally these functions would use a vectorized math library, which PoCL is built with but doesn't use for kernels.

What happens

1. Julia side: SPIRVIntrinsics replaces Base.sin(::Float64) etc. with the OpenCL builtins (@builtin_ccall("sin", …) → _Z3sind), so Julia's own implementation (and its DomainError check) isn't used.
2. PoCL's kernel library: pocl_standalone_jll is built with ENABLE_HOST_CPU_VECTORIZE_LIBMVEC=ON against SLEEF's libsleefgnuabi (Yggdrasil P/pocl/common.jl). In that configuration PoCL's builtins for sin, cos, exp, exp2, log, pow, sqrt, tan are thin wrappers around llvm.sin.f64 and friends, so that "the vectorizer can use vector variants".
3. But the vectorizer never gets them. In PoCL 7.2 (lib/CL/pocl_llvm_wg.cc) the vector library is registered only in initPassManagerForCodeGen, i.e. for code generation after the work-group passes and the loop vectorizer have run. The equivalent TargetLibraryInfo setup in the work-group pipeline (PoCLModulePassManager::build) is under #if 0, and that's still the case on PoCL main.
   • So the loop vectorizer sees llvm.sin.f64 without a vector variant, reports that vectorization is not beneficial, and the kernel ends up with a scalar call to glibc's sin per work-item (nm on the kernel object: U sin, U exp, U log).
   • POCL_DEBUG_LLVM_OPTS=-vector-library=LIBMVEC is ignored ("unknown LLVM option"), so there's no runtime workaround.
4. Not all functions take that path: sincos, sinpi, expm1, log1p aren't in that list and use libclc/SLEEF implementations instead. sinpi gets inlined and is faster than Julia's.

Cost

ns per element, one work-item per element vs. a plain Julia loop over a Vector{Float64}. KA main d048363, pocl_standalone_jll 7.2.2, Julia 1.12, 1 thread, Ryzen 9 5950X (AVX2), pinned; the machine was shared, so ±10%:

@fastmath exp in the kernel takes 2.19 ns, because it uses Julia's exp_impl_fast instead of the builtin.

In Trixi.jl the source terms of elixir_euler_source_terms.jl (sin/cos per node) cost noticeably more on the KA path than on its Polyester path, as part of the gap analysed in #859's investigation.

Options

1. Fix it in PoCL (preferred): register the vector library (createTLII(…, LIBMVEC/SLEEF)) in the work-group pass pipeline as well, so the loop vectorizer can call _ZGVdN4v_sin etc. from SLEEF across work-items. This needs a PoCL patch (upstream, plus a Yggdrasil backport like the other PoCL patches) and should then work for every kernel without changes in Julia. It still only helps where PoCL vectorizes the work-item loop at all, which it currently doesn't for kernels with inner loops or non-inlined calls.
2. Use Julia's implementations on the CPU: override sin/cos/exp/… in the PoCL method table with Base's implementations, which inline and give the vectorizer plain arithmetic. I have a prototype that does this. Float64 results are bitwise identical to Base in 99.98–100% of cases, and PoCL's final kernel then contains no calls.
   • It needs a 64-bit-limb paynehanek, since the SPIR-V backend rejects UInt128.
   • It needs workarounds for OpenCL's ban on recursion, e.g. in ^.
   • Float32 exp loses accuracy (1.12 vs 0.84 ulp), because Base relies on FMAs and contract doesn't survive SPIR-V; that's a separate issue.
3. Bind SLEEF's vector functions directly (as GPU back-ends do with their vendor libraries). That's harder to use from per-work-item code.

With (1) the math stays in PoCL and vectorizes where the loop does. (2) helps kernels that don't vectorize, and fixes the loss of @fastmath. The two aren't exclusive.