```5261 libm should stop using synonyms.h
5298 fabs is 0-sized, confuses dis(1) and others
Reviewed by: Josef 'Jeff' Sipek <jeffpc@josefsipek.net>
Approved by: Gordon Ross <gwr@nexenta.com>
```

 Split Close Expand all Collapse all
```          --- old/usr/src/lib/libm/common/LD/erfl.c
+++ new/usr/src/lib/libm/common/LD/erfl.c
```
 ↓ open down ↓ 24 lines elided ↑ open up ↑
```  25   25  /*
27   27   * Use is subject to license terms.
28   28   */
29   29
30   30  /* long double function erf,erfc (long double x)
31   31   * K.C. Ng, September, 1989.
32   32   *                           x
33   33   *                    2      |\
34   34   *     erf(x)  =  ---------  | exp(-t*t)dt
35      - *                 sqrt(pi) \|
35 + *                 sqrt(pi) \|
36   36   *                           0
37   37   *
38   38   *     erfc(x) =  1-erf(x)
39   39   *
40   40   * method:
41   41   *      Since erf(-x) = -erf(x), we assume x>=0.
42   42   *      For x near 0, we have the expansion
43      - *
43 + *
44   44   *          erf(x) = (2/sqrt(pi))*(x - x^3/3 + x^5/10 - x^7/42 + ....).
45   45   *
46   46   *      Since 2/sqrt(pi) = 1.128379167095512573896158903121545171688,
47   47   *      we use x + x*P(x^2) to approximate erf(x). This formula will
48   48   *      guarantee the error less than one ulp where x is not too far
49   49   *      away from 0. We note that erf(x)=x at x = 0.6174...... After
50      - *      some experiment, we choose the following approximation on
50 + *      some experiment, we choose the following approximation on
51   51   *      interval [0,0.84375].
52   52   *
53   53   *      For x in [0,0.84375]
54   54   *                 2                2        4               40
55   55   *         P =  P(x ) = (p0 + p1 * x + p2 * x + ... + p20 * x  )
56   56   *
57   57   *         erf(x)  = x + x*P
58   58   *         erfc(x) = 1 - erf(x)           if x<=0.25
59   59   *                 = 0.5 + ((0.5-x)-x*P)  if x in [0.25,0.84375]
60   60   *      precision: |P(x^2)-(erf(x)-x)/x| <= 2**-122.50
61   61   *
62      - *      For x in [0.84375,1.25], let s = x - 1, and
62 + *      For x in [0.84375,1.25], let s = x - 1, and
63   63   *      c = 0.84506291151 rounded to single (24 bits)
64   64   *         erf(x)  = c  + P1(s)/Q1(s)
65      - *         erfc(x) = (1-c)  - P1(s)/Q1(s)
65 + *         erfc(x) = (1-c)  - P1(s)/Q1(s)
66   66   *      precision: |P1/Q1 - (erf(x)-c)| <= 2**-118.41
67      - *
68   67   *
69      - *      For x in [1.25,1.75], let s = x - 1.5, and
68 + *
69 + *      For x in [1.25,1.75], let s = x - 1.5, and
70   70   *      c = 0.95478588343 rounded to single (24 bits)
71   71   *         erf(x)  = c  + P2(s)/Q2(s)
72      - *         erfc(x) = (1-c)  - P2(s)/Q2(s)
72 + *         erfc(x) = (1-c)  - P2(s)/Q2(s)
73   73   *      precision: |P1/Q1 - (erf(x)-c)| <= 2**-123.83
74   74   *
75   75   *
76   76   *      For x in [1.75,16/3]
77   77   *         erfc(x) = exp(-x*x)*(1/x)*R1(1/x)/S1(1/x)
78   78   *         erf(x)  = 1 - erfc(x)
79   79   *      precision: absolute error of R1/S1 is bounded by 2**-124.03
80   80   *
81   81   *      For x in [16/3,107]
82   82   *         erfc(x) = exp(-x*x)*(1/x)*R2(1/x)/S2(1/x)
83   83   *         erf(x)  = 1 - erfc(x) (if x>=9 simple return erf(x)=1 with inexact)
84   84   *      precision: absolute error of R2/S2 is bounded by 2**-120.07
85   85   *
86   86   *      Else if inf > x >= 107
87      - *         erf(x)  = 1 with inexact
87 + *         erf(x)  = 1 with inexact
88   88   *         erfc(x) = 0 with underflow
89      - *
89 + *
90   90   *      Special case:
91   91   *         erf(inf)  = 1
92   92   *         erfc(inf) = 0
93   93   */
94   94
95      -#pragma weak erfl = __erfl
96      -#pragma weak erfcl = __erfcl
95 +#pragma weak __erfl = erfl
96 +#pragma weak __erfcl = erfcl
97   97
98   98  #include "libm.h"
99   99  #include "longdouble.h"
100  100
101  101  static long double
102  102  tiny        = 1e-40L,
103  103  nearunfl    = 1e-4000L,
104  104  half        = 0.5L,
105  105  one         = 1.0L,
106  106  onehalf     = 1.5L,
```
 ↓ open down ↓ 18 lines elided ↑ open up ↑
``` 125  125    -6.711366846653439036162105104991433380926e-0012L,
126  126     4.463224090341893165100275380693843116240e-0013L,
127  127    -2.783513452582658245422635662559779162312e-0014L,
128  128     1.634227412586960195251346878863754661546e-0015L,
129  129    -9.060782672889577722765711455623117802795e-0017L,
130  130     4.741341801266246873412159213893613602354e-0018L,
131  131    -2.272417596497826188374846636534317381203e-0019L,
132  132     8.069088733716068462496835658928566920933e-0021L,
133  133  };
134  134
135      -/*
135 +/*
136  136   * Rational erf(x) = ((float)0.84506291151) + P1(x-1)/Q1(x-1) on [0.84375,1.25]
137  137   */
138  138  static long double C1   = (long double)((float)0.84506291151);
139  139  static long double P1[] = {     /*  12 top coeffs */
140  140    -2.362118560752659955654364917390741930316e-0003L,
141  141     4.129623379624420034078926610650759979146e-0001L,
142  142    -3.973857505403547283109417923182669976904e-0002L,
143  143     4.357503184084022439763567513078036755183e-0002L,
144  144     8.015593623388421371247676683754171456950e-0002L,
145  145    -1.034459310403352486685467221776778474602e-0002L,
```
 ↓ open down ↓ 11 lines elided ↑ open up ↑
``` 157  157     1.702332233546316765818144723063881095577e-0001L,
158  158     7.498098377690553934266423088708614219356e-0002L,
159  159     2.050154396918178697056927234366372760310e-0002L,
160  160     7.012988534031999899054782333851905939379e-0003L,
161  161     1.149904787014400354649843451234570731076e-0003L,
162  162     3.185620255011299476196039491205159718620e-0004L,
163  163     1.273405072153008775426376193374105840517e-0005L,
164  164     4.753866999959432971956781228148402971454e-0006L,
165  165    -1.002287602111660026053981728549540200683e-0006L,
166  166  };
167      -/*
167 +/*
168  168   * Rational erf(x) = ((float)0.95478588343) + P2(x-1.5)/Q2(x-1.5)
169  169   * on [1.25,1.75]
170  170   */
171  171  static long double C2   = (long double)((float)0.95478588343);
172  172  static long double P2[] = {     /*  12 top coeffs */
173  173     1.131926304864446730135126164594785863512e-0002L,
174  174     1.273617996967754151544330055186210322832e-0001L,
175  175    -8.169980734667512519897816907190281143423e-0002L,
176  176     9.512267486090321197833634271787944271746e-0002L,
177  177    -2.394251569804872160005274999735914368170e-0002L,
```
 ↓ open down ↓ 89 lines elided ↑ open up ↑
``` 267  267     3.792278350536686111444869752624492443659e+0005L,
268  268     1.257750606950115799965366001773094058720e+0005L,
269  269     3.410830600242369370645608634643620355058e+0004L,
270  270     7.513984469742343134851326863175067271240e+0003L,
271  271     1.313296320593190002554779998138695507840e+0003L,
272  272     1.773972700887629157006326333696896516769e+0002L,
273  273     1.670876451822586800422009013880457094162e+0001L,
274  274     1.000L,
275  275  };
276  276
277      -long double erfl(x)
277 +long double erfl(x)
278  278  long double x;
279  279  {
280  280          long double erfcl(long double),s,y,t;
281  281
282  282          if (!finitel(x)) {
283  283              if (x != x) return x+x;     /* NaN */
284  284              return copysignl(one,x);    /* return +-1.0 is x=Inf */
285  285          }
286  286
287  287          y = fabsl(x);
```
 ↓ open down ↓ 4 lines elided ↑ open up ↑
``` 292  292              return  x+x*t;
293  293          }
294  294          if (y<=1.25L) {
295  295              s = y-one;
296  296              t = C1+__poly_libmq(s,12,P1)/(one+s*__poly_libmq(s,12,Q1));
297  297              return (signbitl(x))? -t: t;
298  298          } else if (y<=1.75L) {
299  299              s = y-onehalf;
300  300              t = C2+__poly_libmq(s,12,P2)/(one+s*__poly_libmq(s,13,Q2));
301  301              return (signbitl(x))? -t: t;
302      -        }
302 +        }
303  303          if (y<=9.0L) t = erfcl(y); else t = tiny;
304  304          return (signbitl(x))? t-one: one-t;
305  305  }
306  306
307  307  long double erfcl(x)
308  308  long double x;
309  309  {
310  310          long double erfl(long double),s,y,t;
311  311
312  312          if (!finitel(x)) {
```
 ↓ open down ↓ 11 lines elided ↑ open up ↑
``` 324  324              return  half+t;
325  325          }
326  326          if (x<=1.25L) {
327  327              s = x-one;
328  328              t = one-C1;
329  329              return t - __poly_libmq(s,12,P1)/(one+s*__poly_libmq(s,12,Q1));
330  330          } else if (x<=1.75L) {
331  331              s = x-onehalf;
332  332              t = one-C2;
333  333              return t - __poly_libmq(s,12,P2)/(one+s*__poly_libmq(s,13,Q2));
334      -        }
334 +        }
335  335          if (x>=107.0L) return nearunfl*nearunfl;                /* underflow */
336  336          else if (x >= L16_3) {
337  337              y = __poly_libmq(x,15,R2);
338  338              t = y/__poly_libmq(x,16,S2);
339  339          } else {
340  340              y = __poly_libmq(x,14,R1);
341  341              t = y/__poly_libmq(x,15,S1);
342  342          }
343  343          /* see comment in ../Q/erfl.c */
344  344          y = x;
345  345          *(int*)&y = 0;
346  346          t *= expl(-y*y)*expl(-(x-y)*(x+y));
347  347          return  t;
348  348  }
```
`XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX`