re PR libfortran/49024 (REAL*16 ERFC_SCALED inaccuracy)

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.

From-SVN: r205151
This commit is contained in:
Francois-Xavier Coudert 2013-11-20 22:18:55 +00:00 committed by François-Xavier Coudert
parent 27b097f8d7
commit 41fd665971
5 changed files with 108 additions and 6 deletions

View File

@ -1,3 +1,8 @@
2013-11-20 Francois-Xavier Coudert <fxcoudert@gcc.gnu.org>
PR libfortran/49024
* gfortran.dg/erf_3.F90: New file.
2013-11-20 Bill Schmidt <wschmidt@linux.vnet.ibm.com>
* gcc.target/powerpc/pr48258-1.c: Skip for little endian.

View File

@ -0,0 +1,44 @@
! { dg-do run }
! { dg-options "-fno-range-check -ffree-line-length-none -O0" }
! { dg-add-options ieee }
!
! Check that simplification functions and runtime library agree on ERF,
! ERFC and ERFC_SCALED, for quadruple-precision.
program test
implicit none
real(kind=16) :: x16
#define CHECK(a) \
x16 = a ; \
call check(erf(real(a,kind=16)), erf(x16)) ; \
call check(erfc(real(a,kind=16)), erfc(x16)) ; \
call check(erfc_scaled(real(a,kind=16)), erfc_scaled(x16))
CHECK(0.0)
CHECK(0.9)
CHECK(1.9)
CHECK(10.)
CHECK(11.)
CHECK(12.)
CHECK(13.)
CHECK(14.)
CHECK(49.)
CHECK(190.)
CHECK(-0.0)
CHECK(-0.9)
CHECK(-1.9)
CHECK(-19.)
CHECK(-190.)
contains
subroutine check (a, b)
real(kind=16), intent(in) :: a, b
print *, abs(a-b) / spacing(a)
if (abs(a - b) > 10 * spacing(a)) call abort
end subroutine
end program test

View File

@ -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

View File

@ -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);
}
}

View File

@ -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?"