// Define which host runtime functions will be used for folding
-// C++ Bessel functions take a floating point as first arguments.
-// Fortran Bessel functions take an integer.
-template<typename HostT> static HostT Bessel_j0(HostT x) {
- return std::cyl_bessel_j(0., x);
-}
-template<typename HostT> static HostT Bessel_y0(HostT x) {
- return std::cyl_neumann(0., x);
-}
-template<typename HostT> static HostT Bessel_j1(HostT x) {
- return std::cyl_bessel_j(1., x);
-}
-template<typename HostT> static HostT Bessel_y1(HostT x) {
- return std::cyl_neumann(1., x);
-}
-template<typename HostT> static HostT Bessel_jn(int n, HostT x) {
- return std::cyl_bessel_j(static_cast<HostT>(n), x);
-}
-template<typename HostT> static HostT Bessel_yn(int n, HostT x) {
- return std::cyl_neumann(static_cast<HostT>(n), x);
-}
+// C++17 defined standard Bessel math functions std::cyl_bessel_j
+// and std::cyl_neumann that can be used for Fortran j and y
+// bessel functions. However, they are not yet implemented in
+// clang libc++ (ok in GNU libstdc++). C maths functions are used
+// in the meantime. They are not C standard but a GNU extension.
+// However, they seem widespread enough to be used.
+
+enum class Bessel { j0, j1, jn, y0, y1, yn };
+template<Bessel, typename T> constexpr auto Sym{0};
+template<> constexpr auto Sym<Bessel::j0, float>{j0f};
+template<> constexpr auto Sym<Bessel::j0, double>{j0};
+template<> constexpr auto Sym<Bessel::j0, long double>{j0l};
+template<> constexpr auto Sym<Bessel::j1, float>{j1f};
+template<> constexpr auto Sym<Bessel::j1, double>{j1};
+template<> constexpr auto Sym<Bessel::j1, long double>{j1l};
+template<> constexpr auto Sym<Bessel::jn, float>{jnf};
+template<> constexpr auto Sym<Bessel::jn, double>{jn};
+template<> constexpr auto Sym<Bessel::jn, long double>{jnl};
+template<> constexpr auto Sym<Bessel::y0, float>{y0f};
+template<> constexpr auto Sym<Bessel::y0, double>{y0};
+template<> constexpr auto Sym<Bessel::y0, long double>{y0l};
+template<> constexpr auto Sym<Bessel::y1, float>{y1f};
+template<> constexpr auto Sym<Bessel::y1, double>{y1};
+template<> constexpr auto Sym<Bessel::y1, long double>{y1l};
+template<> constexpr auto Sym<Bessel::yn, float>{ynf};
+template<> constexpr auto Sym<Bessel::yn, double>{yn};
+template<> constexpr auto Sym<Bessel::yn, long double>{ynl};
template<typename HostT>
static void AddLibmRealHostProcedures(
{"acosh", F{std::acosh}, true}, {"asin", F{std::asin}, true},
{"asinh", F{std::asinh}, true}, {"atan", F{std::atan}, true},
{"atan", F2{std::atan2}, true}, {"atanh", F{std::atanh}, true},
- {"bessel_j0", &Bessel_j0<HostT>, true},
- {"bessel_y0", &Bessel_y0<HostT>, true},
- {"bessel_j1", &Bessel_j1<HostT>, true},
- {"bessel_y1", &Bessel_y1<HostT>, true},
- {"bessel_jn", &Bessel_jn<HostT>, true},
- {"bessel_yn", &Bessel_yn<HostT>, true}, {"cos", F{std::cos}, true},
+ {"bessel_j0", Sym<Bessel::j0, HostT>, true},
+ {"bessel_j1", Sym<Bessel::j1, HostT>, true},
+ {"bessel_jn", Sym<Bessel::jn, HostT>, true},
+ {"bessel_y0", Sym<Bessel::y0, HostT>, true},
+ {"bessel_y1", Sym<Bessel::y1, HostT>, true},
+ {"bessel_yn", Sym<Bessel::yn, HostT>, true}, {"cos", F{std::cos}, true},
{"cosh", F{std::cosh}, true}, {"erf", F{std::erf}, true},
{"erfc", F{std::erfc}, true}, {"exp", F{std::exp}, true},
{"gamma", F{std::tgamma}, true}, {"hypot", F2{std::hypot}, true},
template<L, I, typename> constexpr auto Sym{NoSuchRuntimeSymbol{}};
// Macros to declare fast/relaxed/precise libpgmath variants.
-// Note: std::complex and _complex are layout compatible but only std::comlpex
-// should be used here so that templatized functions work as expected (_Complex
-// and std::complex are different from a type point of view).
#define DECLARE_PGMATH_FAST_REAL(func) \
extern "C" float __fs_##func##_1(float); \
extern "C" double __fd_##func##_1(double); \
template<> constexpr auto Sym<L::F, I::func, double>{__fd_##func##_1};
#define DECLARE_PGMATH_FAST_COMPLEX(func) \
- extern "C" std::complex<float> __fc_##func##_1(std::complex<float>); \
- extern "C" std::complex<double> __fz_##func##_1(std::complex<double>); \
+ extern "C" float _Complex __fc_##func##_1(float _Complex); \
+ extern "C" double _Complex __fz_##func##_1(double _Complex); \
template<> \
constexpr auto Sym<L::F, I::func, std::complex<float>>{__fc_##func##_1}; \
template<> \
template<> constexpr auto Sym<L::P, I::func, double>{__pd_##func##_1};
#define DECLARE_PGMATH_PRECISE_COMPLEX(func) \
- extern "C" std::complex<float> __pc_##func##_1(std::complex<float>); \
- extern "C" std::complex<double> __pz_##func##_1(std::complex<double>); \
+ extern "C" float _Complex __pc_##func##_1(float _Complex); \
+ extern "C" double _Complex __pz_##func##_1(double _Complex); \
template<> \
constexpr auto Sym<L::P, I::func, std::complex<float>>{__pc_##func##_1}; \
template<> \
template<> constexpr auto Sym<L::R, I::func, double>{__rd_##func##_1};
#define DECLARE_PGMATH_RELAXED_COMPLEX(func) \
- extern "C" std::complex<float> __rc_##func##_1(std::complex<float>); \
- extern "C" std::complex<double> __rz_##func##_1(std::complex<double>); \
+ extern "C" float _Complex __rc_##func##_1(float _Complex); \
+ extern "C" double _Complex __rz_##func##_1(double _Complex); \
template<> \
constexpr auto Sym<L::R, I::func, std::complex<float>>{__rc_##func##_1}; \
template<> \
}
}
+// Note: std::complex and _complex are layout compatible but are not guaranteed
+// to be linkage compatible. For instance, on i386, float _Complex is returned
+// by a pair of register but std::complex<float> is returned by structure
+// address. To fix the issue, wrapper around C _Complex functions are defined
+// below.
+template<FuncPointer<float _Complex, float _Complex> func>
+static std::complex<float> ComplexCFuncWrapper(std::complex<float> &arg) {
+ float _Complex res{func(*reinterpret_cast<float _Complex *>(&arg))};
+ return *reinterpret_cast<std::complex<float> *>(&res);
+}
+
+template<FuncPointer<double _Complex, double _Complex> func>
+static std::complex<double> ComplexCFuncWrapper(std::complex<double> &arg) {
+ double _Complex res{func(*reinterpret_cast<double _Complex *>(&arg))};
+ return *reinterpret_cast<std::complex<double> *>(&res);
+}
+
template<L Lib, typename HostT>
static void AddLibpgmathComplexHostProcedures(
HostIntrinsicProceduresLibrary &hostIntrinsicLibrary) {
using CmathF = FuncPointer<CHostT, const CHostT &>;
HostRuntimeIntrinsicProcedure pgmathSymbols[]{
{"abs", FuncPointer<HostT, const CHostT &>{std::abs}, true},
- {"acos", Sym<Lib, I::acos, CHostT>, true},
+ {"acos", ComplexCFuncWrapper<Sym<Lib, I::acos, CHostT>>, true},
{"acosh", CmathF{std::acosh}, true},
- {"asin", Sym<Lib, I::asin, CHostT>, true},
- {"asinh", CmathF{std::asinh}, true}, {"atan", CmathF{std::atan}, true},
+ {"asin", ComplexCFuncWrapper<Sym<Lib, I::asin, CHostT>>, true},
+ {"asinh", CmathF{std::asinh}, true},
+ {"atan", ComplexCFuncWrapper<Sym<Lib, I::atan, CHostT>>, true},
{"atanh", CmathF{std::atanh}, true},
- {"cos", Sym<Lib, I::cos, CHostT>, true},
- {"cosh", Sym<Lib, I::cosh, CHostT>, true},
- {"exp", Sym<Lib, I::exp, CHostT>, true},
- {"log", Sym<Lib, I::log, CHostT>, true},
- {"sin", Sym<Lib, I::sin, CHostT>, true},
- {"sinh", Sym<Lib, I::sinh, CHostT>, true},
- {"sqrt", Sym<Lib, I::sqrt, CHostT>, true},
- {"tan", Sym<Lib, I::tan, CHostT>, true},
- {"tanh", Sym<Lib, I::tanh, CHostT>, true}};
+ {"cos", ComplexCFuncWrapper<Sym<Lib, I::cos, CHostT>>, true},
+ {"cosh", ComplexCFuncWrapper<Sym<Lib, I::cosh, CHostT>>, true},
+ {"exp", ComplexCFuncWrapper<Sym<Lib, I::exp, CHostT>>, true},
+ {"log", ComplexCFuncWrapper<Sym<Lib, I::log, CHostT>>, true},
+ {"sin", ComplexCFuncWrapper<Sym<Lib, I::sin, CHostT>>, true},
+ {"sinh", ComplexCFuncWrapper<Sym<Lib, I::sinh, CHostT>>, true},
+ {"sqrt", ComplexCFuncWrapper<Sym<Lib, I::sqrt, CHostT>>, true},
+ {"tan", ComplexCFuncWrapper<Sym<Lib, I::tan, CHostT>>, true},
+ {"tanh", ComplexCFuncWrapper<Sym<Lib, I::tanh, CHostT>>, true}};
for (auto sym : pgmathSymbols) {
hostIntrinsicLibrary.AddProcedure(std::move(sym));
// available, this needs to be revisited to take it into account. So far,
// default to libpgmath if F18 is built with it.
#if LINK_WITH_LIBPGMATH
- pgmath::InitHostIntrinsicLibraryWithLibpgmath<pgmath::L::P>(*this);
+ // This looks and is stupid for now (until TODO above), but it is needed
+ // to silence clang warnings on unused symbols if all declared pgmath
+ // symbols are not used somewhere.
+ if (true) {
+ pgmath::InitHostIntrinsicLibraryWithLibpgmath<pgmath::L::P>(*this);
+ } else if (false) {
+ pgmath::InitHostIntrinsicLibraryWithLibpgmath<pgmath::L::F>(*this);
+ } else {
+ pgmath::InitHostIntrinsicLibraryWithLibpgmath<pgmath::L::R>(*this);
+ }
#else
InitHostIntrinsicLibraryWithLibm(*this);
#endif