The subject is incorrect, please ignore it, I already fixed it locally. It should be fix bug in erf/erfc.
On Wed, Nov 05, 2014 at 01:14:42PM +0800, Zhigang Gong wrote: > erf/erfc diverge (instead of converging to 1 or 0) for arguments above > about 2. > > This patch is from: > https://bugs.debian.org/cgi-bin/bugreport.cgi?msg=5;filename=Fix-pow-erf-tgamma.patch;att=3;bug=768090 > > Signed-off-by: Zhigang Gong <zhigang.g...@intel.com> > --- > backend/src/ocl_stdlib.tmpl.h | 288 > +++++++++++++++++++++++++++++++++++++++++- > utests/utest_math_gen.py | 16 +-- > 2 files changed, 290 insertions(+), 14 deletions(-) > > diff --git a/backend/src/ocl_stdlib.tmpl.h b/backend/src/ocl_stdlib.tmpl.h > index d6900ee..cacba85 100755 > --- a/backend/src/ocl_stdlib.tmpl.h > +++ b/backend/src/ocl_stdlib.tmpl.h > @@ -2501,12 +2501,6 @@ INLINE_OVERLOADABLE float > __gen_ocl_internal_atan(float x) { > INLINE_OVERLOADABLE float __gen_ocl_internal_atanpi(float x) { > return __gen_ocl_internal_atan(x) / M_PI_F; > } > -INLINE_OVERLOADABLE float __gen_ocl_internal_erf(float x) { > - return M_2_SQRTPI_F * (x - __gen_ocl_pow(x, 3) / 3 + __gen_ocl_pow(x, 5) / > 10 - __gen_ocl_pow(x, 7) / 42 + __gen_ocl_pow(x, 9) / 216); > -} > -INLINE_OVERLOADABLE float __gen_ocl_internal_erfc(float x) { > - return 1 - __gen_ocl_internal_erf(x); > -} > > // XXX work-around PTX profile > #define sqrt native_sqrt > @@ -2726,6 +2720,288 @@ INLINE_OVERLOADABLE float > __gen_ocl_internal_exp(float x) { > return y*twom100; > } > } > + > +/* erf,erfc from glibc s_erff.c -- float version of s_erf.c. > + * Conversion to float by Ian Lance Taylor, Cygnus Support, i...@cygnus.com. > + */ > + > +/* > + * ==================================================== > + * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. > + * > + * Developed at SunPro, a Sun Microsystems, Inc. business. > + * Permission to use, copy, modify, and distribute this > + * software is freely granted, provided that this notice > + * is preserved. > + * ==================================================== > + */ > + > +INLINE_OVERLOADABLE float __gen_ocl_internal_erf(float x) { > +/*...*/ > +const float > +tiny = 1.0e-30, > +half_val= 5.0000000000e-01, /* 0x3F000000 */ > +one = 1.0000000000e+00, /* 0x3F800000 */ > +two = 2.0000000000e+00, /* 0x40000000 */ > + /* c = (subfloat)0.84506291151 */ > +erx = 8.4506291151e-01, /* 0x3f58560b */ > +/* > + * Coefficients for approximation to erf on [0,0.84375] > + */ > +efx = 1.2837916613e-01, /* 0x3e0375d4 */ > +efx8= 1.0270333290e+00, /* 0x3f8375d4 */ > +pp0 = 1.2837916613e-01, /* 0x3e0375d4 */ > +pp1 = -3.2504209876e-01, /* 0xbea66beb */ > +pp2 = -2.8481749818e-02, /* 0xbce9528f */ > +pp3 = -5.7702702470e-03, /* 0xbbbd1489 */ > +pp4 = -2.3763017452e-05, /* 0xb7c756b1 */ > +qq1 = 3.9791721106e-01, /* 0x3ecbbbce */ > +qq2 = 6.5022252500e-02, /* 0x3d852a63 */ > +qq3 = 5.0813062117e-03, /* 0x3ba68116 */ > +qq4 = 1.3249473704e-04, /* 0x390aee49 */ > +qq5 = -3.9602282413e-06, /* 0xb684e21a */ > +/* > + * Coefficients for approximation to erf in [0.84375,1.25] > + */ > +pa0 = -2.3621185683e-03, /* 0xbb1acdc6 */ > +pa1 = 4.1485610604e-01, /* 0x3ed46805 */ > +pa2 = -3.7220788002e-01, /* 0xbebe9208 */ > +pa3 = 3.1834661961e-01, /* 0x3ea2fe54 */ > +pa4 = -1.1089469492e-01, /* 0xbde31cc2 */ > +pa5 = 3.5478305072e-02, /* 0x3d1151b3 */ > +pa6 = -2.1663755178e-03, /* 0xbb0df9c0 */ > +qa1 = 1.0642088205e-01, /* 0x3dd9f331 */ > +qa2 = 5.4039794207e-01, /* 0x3f0a5785 */ > +qa3 = 7.1828655899e-02, /* 0x3d931ae7 */ > +qa4 = 1.2617121637e-01, /* 0x3e013307 */ > +qa5 = 1.3637083583e-02, /* 0x3c5f6e13 */ > +qa6 = 1.1984500103e-02, /* 0x3c445aa3 */ > + /* > + * Coefficients for approximation to erfc in [1.25,1/0.35] > + */ra0 = -9.8649440333e-03, /* 0xbc21a093 */ > +ra1 = -6.9385856390e-01, /* 0xbf31a0b7 */ > +ra2 = -1.0558626175e+01, /* 0xc128f022 */ > +ra3 = -6.2375331879e+01, /* 0xc2798057 */ > +ra4 = -1.6239666748e+02, /* 0xc322658c */ > +ra5 = -1.8460508728e+02, /* 0xc3389ae7 */ > +ra6 = -8.1287437439e+01, /* 0xc2a2932b */ > +ra7 = -9.8143291473e+00, /* 0xc11d077e */ > +sa1 = 1.9651271820e+01, /* 0x419d35ce */ > +sa2 = 1.3765776062e+02, /* 0x4309a863 */ > +sa3 = 4.3456588745e+02, /* 0x43d9486f */ > +sa4 = 6.4538726807e+02, /* 0x442158c9 */ > +sa5 = 4.2900814819e+02, /* 0x43d6810b */ > +sa6 = 1.0863500214e+02, /* 0x42d9451f */ > +sa7 = 6.5702495575e+00, /* 0x40d23f7c */ > +sa8 = -6.0424413532e-02, /* 0xbd777f97 */ > +/* > + * Coefficients for approximation to erfc in [1/.35,28] > + */ > +rb0 = -9.8649431020e-03, /* 0xbc21a092 */ > +rb1 = -7.9928326607e-01, /* 0xbf4c9dd4 */ > +rb2 = -1.7757955551e+01, /* 0xc18e104b */ > +rb3 = -1.6063638306e+02, /* 0xc320a2ea */ > +rb4 = -6.3756646729e+02, /* 0xc41f6441 */ > +rb5 = -1.0250950928e+03, /* 0xc480230b */ > +rb6 = -4.8351919556e+02, /* 0xc3f1c275 */ > +sb1 = 3.0338060379e+01, /* 0x41f2b459 */ > +sb2 = 3.2579251099e+02, /* 0x43a2e571 */ > +sb3 = 1.5367296143e+03, /* 0x44c01759 */ > +sb4 = 3.1998581543e+03, /* 0x4547fdbb */ > +sb5 = 2.5530502930e+03, /* 0x451f90ce */ > +sb6 = 4.7452853394e+02, /* 0x43ed43a7 */ > +sb7 = -2.2440952301e+01; /* 0xc1b38712 */ > + > + int hx,ix,i; > + float R,S,P,Q,s,y,z,r; > + GEN_OCL_GET_FLOAT_WORD(hx,x); > + ix = hx&0x7fffffff; > + if(ix>=0x7f800000) { /* erf(nan)=nan */ > + i = ((unsigned int)hx>>31)<<1; > + return (float)(1-i)+one/x; /* erf(+-inf)=+-1 */ > + } > + > + if(ix < 0x3f580000) { /* |x|<0.84375 */ > + if(ix < 0x31800000) { /* |x|<2**-28 */ > + if (ix < 0x04000000) > + /*avoid underflow */ > + return (float)0.125*((float)8.0*x+efx8*x); > + return x + efx*x; > + } > + z = x*x; > + r = pp0+z*(pp1+z*(pp2+z*(pp3+z*pp4))); > + s = one+z*(qq1+z*(qq2+z*(qq3+z*(qq4+z*qq5)))); > + y = r/s; > + return x + x*y; > + } > + if(ix < 0x3fa00000) { /* 0.84375 <= |x| < 1.25 */ > + s = __gen_ocl_internal_fabs(x)-one; > + P = pa0+s*(pa1+s*(pa2+s*(pa3+s*(pa4+s*(pa5+s*pa6))))); > + Q = one+s*(qa1+s*(qa2+s*(qa3+s*(qa4+s*(qa5+s*qa6))))); > + if(hx>=0) return erx + P/Q; else return -erx - P/Q; > + } > + if (ix >= 0x40c00000) { /* inf>|x|>=6 */ > + if(hx>=0) return one-tiny; else return tiny-one; > + } > + x = __gen_ocl_internal_fabs(x); > + s = one/(x*x); > + if(ix< 0x4036DB6E) { /* |x| < 1/0.35 */ > + R=ra0+s*(ra1+s*(ra2+s*(ra3+s*(ra4+s*( > + ra5+s*(ra6+s*ra7)))))); > + S=one+s*(sa1+s*(sa2+s*(sa3+s*(sa4+s*( > + sa5+s*(sa6+s*(sa7+s*sa8))))))); > + } else { /* |x| >= 1/0.35 */ > + R=rb0+s*(rb1+s*(rb2+s*(rb3+s*(rb4+s*( > + rb5+s*rb6))))); > + S=one+s*(sb1+s*(sb2+s*(sb3+s*(sb4+s*( > + sb5+s*(sb6+s*sb7)))))); > + } > + GEN_OCL_GET_FLOAT_WORD(ix,x); > + GEN_OCL_SET_FLOAT_WORD(z,ix&0xfffff000); > + r = > __gen_ocl_internal_exp(-z*z-(float)0.5625)*__gen_ocl_internal_exp((z-x)*(z+x)+R/S); > + if(hx>=0) return one-r/x; else return r/x-one; > +} > +INLINE_OVERLOADABLE float __gen_ocl_internal_erfc(float x) { > +/*...*/ > +const float > +tiny = 1.0e-30, > +half_val= 5.0000000000e-01, /* 0x3F000000 */ > +one = 1.0000000000e+00, /* 0x3F800000 */ > +two = 2.0000000000e+00, /* 0x40000000 */ > + /* c = (subfloat)0.84506291151 */ > +erx = 8.4506291151e-01, /* 0x3f58560b */ > +/* > + * Coefficients for approximation to erf on [0,0.84375] > + */ > +efx = 1.2837916613e-01, /* 0x3e0375d4 */ > +efx8= 1.0270333290e+00, /* 0x3f8375d4 */ > +pp0 = 1.2837916613e-01, /* 0x3e0375d4 */ > +pp1 = -3.2504209876e-01, /* 0xbea66beb */ > +pp2 = -2.8481749818e-02, /* 0xbce9528f */ > +pp3 = -5.7702702470e-03, /* 0xbbbd1489 */ > +pp4 = -2.3763017452e-05, /* 0xb7c756b1 */ > +qq1 = 3.9791721106e-01, /* 0x3ecbbbce */ > +qq2 = 6.5022252500e-02, /* 0x3d852a63 */ > +qq3 = 5.0813062117e-03, /* 0x3ba68116 */ > +qq4 = 1.3249473704e-04, /* 0x390aee49 */ > +qq5 = -3.9602282413e-06, /* 0xb684e21a */ > +/* > + * Coefficients for approximation to erf in [0.84375,1.25] > + */ > +pa0 = -2.3621185683e-03, /* 0xbb1acdc6 */ > +pa1 = 4.1485610604e-01, /* 0x3ed46805 */ > +pa2 = -3.7220788002e-01, /* 0xbebe9208 */ > +pa3 = 3.1834661961e-01, /* 0x3ea2fe54 */ > +pa4 = -1.1089469492e-01, /* 0xbde31cc2 */ > +pa5 = 3.5478305072e-02, /* 0x3d1151b3 */ > +pa6 = -2.1663755178e-03, /* 0xbb0df9c0 */ > +qa1 = 1.0642088205e-01, /* 0x3dd9f331 */ > +qa2 = 5.4039794207e-01, /* 0x3f0a5785 */ > +qa3 = 7.1828655899e-02, /* 0x3d931ae7 */ > +qa4 = 1.2617121637e-01, /* 0x3e013307 */ > +qa5 = 1.3637083583e-02, /* 0x3c5f6e13 */ > +qa6 = 1.1984500103e-02, /* 0x3c445aa3 */ > + /* > + * Coefficients for approximation to erfc in [1.25,1/0.35] > + */ra0 = -9.8649440333e-03, /* 0xbc21a093 */ > +ra1 = -6.9385856390e-01, /* 0xbf31a0b7 */ > +ra2 = -1.0558626175e+01, /* 0xc128f022 */ > +ra3 = -6.2375331879e+01, /* 0xc2798057 */ > +ra4 = -1.6239666748e+02, /* 0xc322658c */ > +ra5 = -1.8460508728e+02, /* 0xc3389ae7 */ > +ra6 = -8.1287437439e+01, /* 0xc2a2932b */ > +ra7 = -9.8143291473e+00, /* 0xc11d077e */ > +sa1 = 1.9651271820e+01, /* 0x419d35ce */ > +sa2 = 1.3765776062e+02, /* 0x4309a863 */ > +sa3 = 4.3456588745e+02, /* 0x43d9486f */ > +sa4 = 6.4538726807e+02, /* 0x442158c9 */ > +sa5 = 4.2900814819e+02, /* 0x43d6810b */ > +sa6 = 1.0863500214e+02, /* 0x42d9451f */ > +sa7 = 6.5702495575e+00, /* 0x40d23f7c */ > +sa8 = -6.0424413532e-02, /* 0xbd777f97 */ > +/* > + * Coefficients for approximation to erfc in [1/.35,28] > + */ > +rb0 = -9.8649431020e-03, /* 0xbc21a092 */ > +rb1 = -7.9928326607e-01, /* 0xbf4c9dd4 */ > +rb2 = -1.7757955551e+01, /* 0xc18e104b */ > +rb3 = -1.6063638306e+02, /* 0xc320a2ea */ > +rb4 = -6.3756646729e+02, /* 0xc41f6441 */ > +rb5 = -1.0250950928e+03, /* 0xc480230b */ > +rb6 = -4.8351919556e+02, /* 0xc3f1c275 */ > +sb1 = 3.0338060379e+01, /* 0x41f2b459 */ > +sb2 = 3.2579251099e+02, /* 0x43a2e571 */ > +sb3 = 1.5367296143e+03, /* 0x44c01759 */ > +sb4 = 3.1998581543e+03, /* 0x4547fdbb */ > +sb5 = 2.5530502930e+03, /* 0x451f90ce */ > +sb6 = 4.7452853394e+02, /* 0x43ed43a7 */ > +sb7 = -2.2440952301e+01; /* 0xc1b38712 */ > + int hx,ix; > + float R,S,P,Q,s,y,z,r; > + GEN_OCL_GET_FLOAT_WORD(hx,x); > + ix = hx&0x7fffffff; > + if(ix>=0x7f800000) { /* erfc(nan)=nan */ > + /* erfc(+-inf)=0,2 */ > + return (float)(((unsigned int)hx>>31)<<1)+one/x; > + } > + > + if(ix < 0x3f580000) { /* |x|<0.84375 */ > + if(ix < 0x23800000) /* |x|<2**-56 */ > + return one-x; > + z = x*x; > + r = pp0+z*(pp1+z*(pp2+z*(pp3+z*pp4))); > + s = one+z*(qq1+z*(qq2+z*(qq3+z*(qq4+z*qq5)))); > + y = r/s; > + if(hx < 0x3e800000) { /* x<1/4 */ > + return one-(x+x*y); > + } else { > + r = x*y; > + r += (x-half_val); > + return half_val - r ; > + } > + } > + if(ix < 0x3fa00000) { /* 0.84375 <= |x| < 1.25 */ > + s = __gen_ocl_internal_fabs(x)-one; > + P = pa0+s*(pa1+s*(pa2+s*(pa3+s*(pa4+s*(pa5+s*pa6))))); > + Q = one+s*(qa1+s*(qa2+s*(qa3+s*(qa4+s*(qa5+s*qa6))))); > + if(hx>=0) { > + z = one-erx; return z - P/Q; > + } else { > + z = erx+P/Q; return one+z; > + } > + } > + if (ix < 0x41e00000) { /* |x|<28 */ > + x = __gen_ocl_internal_fabs(x); > + s = one/(x*x); > + if(ix< 0x4036DB6D) { /* |x| < 1/.35 ~ 2.857143*/ > + R=ra0+s*(ra1+s*(ra2+s*(ra3+s*(ra4+s*( > + ra5+s*(ra6+s*ra7)))))); > + S=one+s*(sa1+s*(sa2+s*(sa3+s*(sa4+s*( > + sa5+s*(sa6+s*(sa7+s*sa8))))))); > + } else { /* |x| >= 1/.35 ~ 2.857143 */ > + if(hx<0&&ix>=0x40c00000) return two-tiny;/* x < -6 */ > + R=rb0+s*(rb1+s*(rb2+s*(rb3+s*(rb4+s*( > + rb5+s*rb6))))); > + S=one+s*(sb1+s*(sb2+s*(sb3+s*(sb4+s*( > + sb5+s*(sb6+s*sb7)))))); > + } > + GEN_OCL_GET_FLOAT_WORD(ix,x); > + GEN_OCL_SET_FLOAT_WORD(z,ix&0xffffe000); > + r = __gen_ocl_internal_exp(-z*z-(float)0.5625)* > + __gen_ocl_internal_exp((z-x)*(z+x)+R/S); > + if(hx>0) { > + float ret = r/x; > + return ret; > + } else > + return two-r/x; > + } else { > + if(hx>0) { > + return tiny*tiny; > + } else > + return two-tiny; > + } > +} > + > INLINE_OVERLOADABLE float __gen_ocl_internal_fmod (float x, float y) { > //return x-y*__gen_ocl_rndz(x/y); > float one = 1.0; > diff --git a/utests/utest_math_gen.py b/utests/utest_math_gen.py > index 30a9b24..69b42a1 100755 > --- a/utests/utest_math_gen.py > +++ b/utests/utest_math_gen.py > @@ -217,16 +217,16 @@ static float atanpi(float x){ > cospiUtests = > func('cospi','cospi',[cospi_input_type],cospi_output_type,[cospi_input_values],'2 > * FLT_ULP',cospi_cpu_func) > > # ##### gentype erf(gentype) > -# erf_input_values = base_input_values > -# erf_input_type = ['float','float2','float4','float8','float16'] > -# erf_output_type = ['float','float2','float4','float8','float16'] > -# erfUtests = > func('erf','erf',[erf_input_type],erf_output_type,[erf_input_values],'16 * > FLT_ULP') > + erf_input_values = base_input_values > + erf_input_type = ['float','float2','float4','float8','float16'] > + erf_output_type = ['float','float2','float4','float8','float16'] > + erfUtests = > func('erf','erf',[erf_input_type],erf_output_type,[erf_input_values],'16 * > FLT_ULP') > > # ##### gentype erfc(gentype) > -# erfc_input_values = base_input_values > -# erfc_input_type = ['float','float2','float4','float8','float16'] > -# erfc_output_type = ['float','float2','float4','float8','float16'] > -# erfcUtests = > func('erfc','erfc',[erfc_input_type],erfc_output_type,[erfc_input_values],'16 > * FLT_ULP') > + erfc_input_values = base_input_values > + erfc_input_type = ['float','float2','float4','float8','float16'] > + erfc_output_type = ['float','float2','float4','float8','float16'] > + erfcUtests = > func('erfc','erfc',[erfc_input_type],erfc_output_type,[erfc_input_values],'16 > * FLT_ULP') > > ##### gentype exp(gentype x) > exp_input_values = base_input_values > -- > 1.8.3.2 > > _______________________________________________ > Beignet mailing list > Beignet@lists.freedesktop.org > http://lists.freedesktop.org/mailman/listinfo/beignet _______________________________________________ Beignet mailing list Beignet@lists.freedesktop.org http://lists.freedesktop.org/mailman/listinfo/beignet