mirror of
https://github.com/bolero-MURAKAMI/Sprout
synced 2024-11-12 21:09:01 +00:00
177 lines
7.6 KiB
C++
177 lines
7.6 KiB
C++
/*=============================================================================
|
|
Copyright (c) 2011-2019 Bolero MURAKAMI
|
|
https://github.com/bolero-MURAKAMI/Sprout
|
|
|
|
Distributed under the Boost Software License, Version 1.0. (See accompanying
|
|
file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
|
|
=============================================================================*/
|
|
#ifndef SPROUT_MATH_ERF_HPP
|
|
#define SPROUT_MATH_ERF_HPP
|
|
|
|
#include <type_traits>
|
|
#include <sprout/config.hpp>
|
|
#include <sprout/limits.hpp>
|
|
#include <sprout/detail/pow.hpp>
|
|
#include <sprout/math/detail/config.hpp>
|
|
#include <sprout/math/detail/float_compute.hpp>
|
|
#include <sprout/math/isnan.hpp>
|
|
#include <sprout/math/itrunc.hpp>
|
|
#include <sprout/math/fractional_part.hpp>
|
|
#include <sprout/type_traits/enabler_if.hpp>
|
|
|
|
namespace sprout {
|
|
namespace math {
|
|
namespace detail {
|
|
#if SPROUT_USE_BUILTIN_CMATH_FUNCTION
|
|
inline SPROUT_CONSTEXPR float
|
|
builtin_erf(float x) {
|
|
return __builtin_erff(x);
|
|
}
|
|
inline SPROUT_CONSTEXPR double
|
|
builtin_erf(double x) {
|
|
return __builtin_erf(x);
|
|
}
|
|
inline SPROUT_CONSTEXPR long double
|
|
builtin_erf(long double x) {
|
|
return __builtin_erfl(x);
|
|
}
|
|
#endif
|
|
|
|
template<typename T>
|
|
inline SPROUT_CONSTEXPR T
|
|
erf_impl_3(T x, T y) {
|
|
return x < 0 ? -y
|
|
: y
|
|
;
|
|
}
|
|
template<typename T>
|
|
inline SPROUT_CONSTEXPR T
|
|
erf_impl_2_b(T x, T t, int k) {
|
|
return sprout::math::detail::erf_impl_3(
|
|
x,
|
|
T(1) - sprout::detail::pow_n(
|
|
k == 0 ? (((((((((((
|
|
-2.9734388465e-10 * t + 2.69776334046e-9) * t
|
|
+ -6.40788827665e-9) * t + -1.6678201321e-8) * t + -2.1854388148686e-7) * t
|
|
+ 2.66246030457984e-6) * t + 1.612722157047886e-5) * t + -2.5616361025506629e-4) * t
|
|
+ 1.5380842432375365e-4) * t + 0.00815533022524927908) * t + -0.01402283663896319337) * t
|
|
+ -0.19746892495383021487) * t + 0.71511720328842845913
|
|
: k == 1 ? (((((((((((
|
|
-1.951073787e-11 * t + -3.2302692214e-10) * t
|
|
+ 5.22461866919e-9) * t + 3.42940918551e-9) * t + -3.5772874310272e-7) * t
|
|
+ 1.9999935792654e-7) * t + 2.687044575042908e-5) * t + -1.1843240273775776e-4) * t
|
|
+ -8.0991728956032271e-4) * t + 0.00661062970502241174) * t + 0.00909530922354827295) * t
|
|
+ -0.2016007277849101314) * t + 0.51169696718727644908
|
|
: k == 2 ? (((((((((((
|
|
3.147682272e-11 * t + -4.8465972408e-10) * t
|
|
+ 6.3675740242e-10) * t + 3.377623323271e-8) * t + -1.5451139637086e-7) * t
|
|
+ -2.03340624738438e-6) * t + 1.947204525295057e-5) * t + 2.854147231653228e-5) * t
|
|
+ -0.00101565063152200272) * t + 0.00271187003520095655) * t + 0.02328095035422810727) * t
|
|
+ -0.16725021123116877197) * t + 0.32490054966649436974
|
|
: k == 3 ? (((((((((((
|
|
2.31936337e-11 * t + -6.303206648e-11) * t
|
|
+ -2.64888267434e-9) * t + 2.050708040581e-8) * t + 1.1371857327578e-7) * t
|
|
+ -2.11211337219663e-6) * t + 3.68797328322935e-6) * t + 9.823686253424796e-5) * t
|
|
+ -6.5860243990455368e-4) * t + -7.5285814895230877e-4) * t + 0.02585434424202960464) * t
|
|
+ -0.11637092784486193258) * t + 0.18267336775296612024
|
|
: (((((((((((
|
|
-3.67789363e-12 * t + 2.0876046746e-10) * t
|
|
+ -1.93319027226e-9) * t + -4.35953392472e-9) * t + 1.8006992266137e-7) * t
|
|
+ -7.8441223763969e-7) * t + -6.75407647949153e-6) * t + 8.428418334440096e-5) * t
|
|
+ -1.7604388937031815e-4) * t + -0.0023972961143507161) * t + 0.0206412902387602297) * t
|
|
+ -0.06905562880005864105) * t + 0.09084526782065478489
|
|
,
|
|
16
|
|
)
|
|
);
|
|
}
|
|
template<typename T>
|
|
inline SPROUT_CONSTEXPR T
|
|
erf_impl_2_a_1(T x, T w, T t, int k) {
|
|
return sprout::math::detail::erf_impl_3(
|
|
x,
|
|
k == 0 ? ((((((((((((
|
|
5.958930743e-11 * t + -1.13739022964e-9) * t
|
|
+ 1.466005199839e-8) * t + -1.635035446196e-7) * t + 1.6461004480962e-6) * t
|
|
+ -1.492559551950604e-5) * t + 1.2055331122299265e-4) * t + -8.548326981129666e-4) * t
|
|
+ 0.00522397762482322257) * t + -0.0268661706450773342) * t + 0.11283791670954881569) * t
|
|
+ -0.37612638903183748117) * t + 1.12837916709551257377) * w
|
|
: k == 1 ? ((((((((((((
|
|
2.372510631e-11 * t + -4.5493253732e-10) * t
|
|
+ 5.90362766598e-9) * t + -6.642090827576e-8) * t + 6.7595634268133e-7) * t
|
|
+ -6.21188515924e-6) * t + 5.10388300970969e-5) * t + -3.7015410692956173e-4) * t
|
|
+ 0.00233307631218880978) * t + -0.0125498847718219221) * t + 0.05657061146827041994) * t
|
|
+ -0.2137966477645600658) * t + 0.84270079294971486929) * w
|
|
: k == 2 ? ((((((((((((
|
|
9.49905026e-12 * t + -1.8310229805e-10) * t
|
|
+ 2.39463074e-9) * t + -2.721444369609e-8) * t + 2.8045522331686e-7) * t
|
|
+ -2.61830022482897e-6) * t + 2.195455056768781e-5) * t + -1.6358986921372656e-4) * t
|
|
+ 0.00107052153564110318) * t + -0.00608284718113590151) * t + 0.02986978465246258244) * t
|
|
+ -0.13055593046562267625) * t + 0.67493323603965504676) * w
|
|
: k == 3 ? ((((((((((((
|
|
3.82722073e-12 * t + -7.421598602e-11) * t
|
|
+ 9.793057408e-10) * t + -1.126008898854e-8) * t + 1.1775134830784e-7) * t
|
|
+ -1.1199275838265e-6) * t + 9.62023443095201e-6) * t + -7.404402135070773e-5) * t
|
|
+ 5.0689993654144881e-4) * t + -0.00307553051439272889) * t + 0.01668977892553165586) * t
|
|
+ -0.08548534594781312114) * t + 0.56909076642393639985) * w
|
|
: ((((((((((((
|
|
1.55296588e-12 * t + -3.032205868e-11) * t
|
|
+ 4.0424830707e-10) * t + -4.71135111493e-9) * t + 5.011915876293e-8) * t
|
|
+ -4.8722516178974e-7) * t + 4.30683284629395e-6) * t + -3.445026145385764e-5) * t
|
|
+ 2.4879276133931664e-4) * t + -0.00162940941748079288) * t + 0.00988786373932350462) * t
|
|
+ -0.05962426839442303805) * t + 0.49766113250947636708) * w
|
|
);
|
|
}
|
|
template<typename T>
|
|
inline SPROUT_CONSTEXPR T
|
|
erf_impl_2_a(T x, T w, T t) {
|
|
return sprout::math::detail::erf_impl_2_a_1(x, w, sprout::math::fractional_part(t), sprout::math::itrunc(t));
|
|
}
|
|
template<typename T>
|
|
inline SPROUT_CONSTEXPR T
|
|
erf_impl_1(T x, T w) {
|
|
return w < T(2.2) ? sprout::math::detail::erf_impl_2_a(x, w, w * w)
|
|
: w < T(6.9) ? sprout::math::detail::erf_impl_2_b(x, sprout::math::fractional_part(w), sprout::math::itrunc(w) - 2)
|
|
: sprout::math::detail::erf_impl_3(x, T(1))
|
|
;
|
|
}
|
|
template<typename T>
|
|
inline SPROUT_CONSTEXPR T
|
|
erf_impl(T x) {
|
|
return sprout::math::detail::erf_impl_1(x, x < 0 ? -x : x);
|
|
}
|
|
} // namespace detail
|
|
//
|
|
// erf
|
|
//
|
|
template<
|
|
typename FloatType,
|
|
typename sprout::enabler_if<std::is_floating_point<FloatType>::value>::type = sprout::enabler
|
|
>
|
|
inline SPROUT_CONSTEXPR FloatType
|
|
erf(FloatType x) {
|
|
return sprout::math::isnan(x) ? x
|
|
: x == sprout::numeric_limits<FloatType>::infinity() ? FloatType(1)
|
|
: x == -sprout::numeric_limits<FloatType>::infinity() ? FloatType(-1)
|
|
#if SPROUT_USE_BUILTIN_CMATH_FUNCTION
|
|
: sprout::math::detail::builtin_erf(x)
|
|
#else
|
|
: x == 0 ? x
|
|
: static_cast<FloatType>(sprout::math::detail::erf_impl(static_cast<typename sprout::math::detail::float_compute<FloatType>::type>(x)))
|
|
#endif
|
|
;
|
|
}
|
|
template<
|
|
typename IntType,
|
|
typename sprout::enabler_if<std::is_integral<IntType>::value>::type = sprout::enabler
|
|
>
|
|
inline SPROUT_CONSTEXPR double
|
|
erf(IntType x) {
|
|
return sprout::math::erf(static_cast<double>(x));
|
|
}
|
|
} // namespace math
|
|
|
|
using sprout::math::erf;
|
|
} // namespace sprout
|
|
|
|
#endif // #ifndef SPROUT_MATH_ERF_HPP
|