diff options
| author | fxcoudert <fxcoudert@138bc75d-0d04-0410-961f-82ee72b054a4> | 2013-11-20 22:18:55 +0000 |
|---|---|---|
| committer | fxcoudert <fxcoudert@138bc75d-0d04-0410-961f-82ee72b054a4> | 2013-11-20 22:18:55 +0000 |
| commit | 62776fa2f34fc239ad60e5c99ea53271dca2045e (patch) | |
| tree | 9e450558efefb68f1ce323672b0891b1418ada2e /libgfortran | |
| parent | 637a491d84b1891636d04b2dec4f34c6b72a28ce (diff) | |
| download | gcc-62776fa2f34fc239ad60e5c99ea53271dca2045e.tar.gz | |
PR libfortran/49024
* intrinsics/erfc_scaled.c (erfc_scaled_r16): New function.
* intrinsics/erfc_scaled_inc.c: Do not provide quadruple
precision variant.
* gfortran.dg/erf_3.F90: New file.
git-svn-id: svn+ssh://gcc.gnu.org/svn/gcc/trunk@205151 138bc75d-0d04-0410-961f-82ee72b054a4
Diffstat (limited to 'libgfortran')
| -rw-r--r-- | libgfortran/ChangeLog | 7 | ||||
| -rw-r--r-- | libgfortran/intrinsics/erfc_scaled.c | 53 | ||||
| -rw-r--r-- | libgfortran/intrinsics/erfc_scaled_inc.c | 5 |
3 files changed, 59 insertions, 6 deletions
diff --git a/libgfortran/ChangeLog b/libgfortran/ChangeLog index fcbc5483472..eeeaa04b3c3 100644 --- a/libgfortran/ChangeLog +++ b/libgfortran/ChangeLog @@ -1,3 +1,10 @@ +2013-11-20 Francois-Xavier Coudert <fxcoudert@gcc.gnu.org> + + PR libfortran/49024 + * intrinsics/erfc_scaled.c (erfc_scaled_r16): New function. + * intrinsics/erfc_scaled_inc.c: Do not provide quadruple + precision variant. + 2013-11-18 Francois-Xavier Coudert <fxcoudert@gcc.gnu.org> PR libfortran/51828 diff --git a/libgfortran/intrinsics/erfc_scaled.c b/libgfortran/intrinsics/erfc_scaled.c index 1c58a08938e..1f8c778eb65 100644 --- a/libgfortran/intrinsics/erfc_scaled.c +++ b/libgfortran/intrinsics/erfc_scaled.c @@ -45,8 +45,59 @@ see the files COPYING3 and COPYING.RUNTIME respectively. If not, see #include "erfc_scaled_inc.c" #endif -#ifdef HAVE_GFC_REAL_16 +#if defined(HAVE_GFC_REAL_16) && defined(GFC_REAL_16_IS_LONG_DOUBLE) #undef KIND #define KIND 16 #include "erfc_scaled_inc.c" #endif + + +/* For quadruple-precision (__float128), netlib's implementation is + not accurate enough. We provide another one. */ + + +extern GFC_REAL_16 erfc_scaled_r16 (GFC_REAL_16); +export_proto(erfc_scaled_r16); + + +GFC_REAL_16 +erfc_scaled_r16 (GFC_REAL_16 x) +{ + if (x < -106.566990228185312813205074546585730Q) + { + return __builtin_infq(); + } + if (x < 12) + { + /* Compute directly as ERFC_SCALED(x) = ERFC(x) * EXP(X**2). + This is not perfect, but much better than netlib. */ + return erfcq(x) * expq(x * x); + } + else + { + /* Calculate ERFC_SCALED(x) using a power series in 1/x: + ERFC_SCALED(x) = 1 / (x * sqrt(pi)) + * (1 + Sum_n (-1)**n * (1 * 3 * 5 * ... * (2n-1)) + / (2 * x**2)**n) + */ + GFC_REAL_16 sum = 0, oldsum; + GFC_REAL_16 inv2x2 = 1 / (2 * x * x); + GFC_REAL_16 fac = 1; + int n = 1; + + while (n < 200) + { + fac *= - (2*n - 1) * inv2x2; + oldsum = sum; + sum += fac; + + if (sum == oldsum) + break; + + n++; + } + + return (1 + sum) / x * (M_2_SQRTPIq / 2); + } +} + diff --git a/libgfortran/intrinsics/erfc_scaled_inc.c b/libgfortran/intrinsics/erfc_scaled_inc.c index 57a6b71f995..107d91a6c9d 100644 --- a/libgfortran/intrinsics/erfc_scaled_inc.c +++ b/libgfortran/intrinsics/erfc_scaled_inc.c @@ -48,11 +48,6 @@ see the files COPYING3 and COPYING.RUNTIME respectively. If not, see # define TRUNC(x) truncl(x) # endif -#elif (KIND == 16 && defined(GFC_REAL_16_IS_FLOAT128)) - -# define EXP(x) expq(x) -# define TRUNC(x) truncq(x) - #else # error "What exactly is it that you want me to do?" |