cyl_bessel_j functions (on device)

[ax3l]

In ECP WarpX, we would like to use the cylindrical Bessel function of the first kind on device.

They are defined in the C++17 standard: https://en.cppreference.com/w/cpp/numeric/special_functions/cyl_bessel_j

But compiling with those device functions errors as of the latest oneAPI/DPC++ compiler out with:

/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/tr1/bessel_function.tcc:495:14: error: SYCL kernel cannot call an undefined function without SYCL_EXTERNAL attribute
        std::__throw_domain_error(__N("Bad argument "
             ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/bits/functexcept.h:66:3: note: '__throw_domain_error' declared here
  __throw_domain_error(const char*) __attribute__((__noreturn__));
  ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/tr1/bessel_function.tcc:492:5: note: called by '__cyl_bessel_j<float>'
    __cyl_bessel_j(_Tp __nu, _Tp __x)
    ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/tr1/bessel_function.tcc:207:14: error: SYCL kernel cannot call an undefined function without SYCL_EXTERNAL attribute
        std::__throw_runtime_error(__N("Argument x too large in __bessel_jn; "
             ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/bits/functexcept.h:82:3: note: '__throw_runtime_error' declared here
  __throw_runtime_error(const char*) __attribute__((__noreturn__));
  ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/tr1/bessel_function.tcc:138:5: note: called by '__bessel_jn<float>'
    __bessel_jn(_Tp __nu, _Tp __x,
    ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/tr1/bessel_function.tcc:264:18: error: SYCL kernel cannot call an undefined function without SYCL_EXTERNAL attribute
            std::__throw_runtime_error(__N("Bessel y series failed to converge "
                 ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/bits/functexcept.h:82:3: note: '__throw_runtime_error' declared here
  __throw_runtime_error(const char*) __attribute__((__noreturn__));
  ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/tr1/bessel_function.tcc:138:5: note: called by '__bessel_jn<float>'
    __bessel_jn(_Tp __nu, _Tp __x,
    ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/tr1/bessel_function.tcc:315:18: error: SYCL kernel cannot call an undefined function without SYCL_EXTERNAL attribute
            std::__throw_runtime_error(__N("Lentz's method failed "
                 ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/bits/functexcept.h:82:3: note: '__throw_runtime_error' declared here
  __throw_runtime_error(const char*) __attribute__((__noreturn__));
  ^
/usr/lib/gcc/x86_64-linux-gnu/10/../../../../include/c++/10/tr1/bessel_function.tcc:138:5: note: called by '__bessel_jn<float>'
    __bessel_jn(_Tp __nu, _Tp __x,

It seems there is a similar error for jn/jnf as defined here: https://software.intel.com/content/www/us/en/develop/documentation/oneapi-dpcpp-cpp-compiler-dev-guide-and-reference/top/optimization-and-programming-guide/intel-math-library/math-functions/special-functions.html

I cannot spot jn(f) in the sycl:: namespace here: https://github.com/intel/llvm/blob/sycl/sycl/doc/extensions/C-CXX-StandardLibrary/C-CXX-StandardLibrary.rst

$ dpcpp --version
Intel(R) oneAPI DPC++ Compiler 2021.1.2 (2020.10.0.1214)
Target: x86_64-unknown-linux-gnu
Thread model: posix
InstalledDir: /opt/intel/oneapi/compiler/2021.1.2/linux/bin

[rscohn2]

@mkinsner : Is this an issue with the spec or the implementation?

[ax3l]

I think the request for std::cyl_bessel_j(f|l) might be a spec feature request and the jn(f) issue might be an implementation issue. Please don't hesitate to redirect :-)

[mkinsner]

oneDPL lists some functions as callable from device code, documented at:

Tested Standard C++ APIs

@akukanov may have more info on plans for bessel variants.

jn/jnf is probably an impl bug. @bader, thoughts?

[bader]

jn/jnf is probably an impl bug. @bader, thoughts?

I haven't seen requests to support this functions before, so I'm not sure if we can classify this as an implementation bug. So far, the math functions enabled in DPC++ have been natively supported by low-level runtimes (e.g. OpenCL), so usually it doesn't take much time to write a DPC++ wrapper for an OpenCL built-in function. This case seems to be different, jn/jnf aren't standard OpenCL built-ins.

I suppose the functions defined here: https://software.intel.com/content/www/us/en/develop/documentation/oneapi-dpcpp-cpp-compiler-dev-guide-and-reference/top/optimization-and-programming-guide/intel-math-library/math-functions/special-functions.html are supported on CPU only with disabled SYCL mode.