Author: kargl
Date: Mon Sep 15 23:21:57 2014
New Revision: 271651
URL: http://svnweb.freebsd.org/changeset/base/271651
Log:
* Makefile:
. Hook e_lgammal[_r].c to the build.
. Create man page links for lgammal[-r].3.
* Symbol.map:
. Sort lgammal to its rightful place.
. Add FBSD_1.4 section for the new lgamal_r symbol.
* ld128/e_lgammal_r.c:
. 128-bit implementataion of lgammal_r().
* ld80/e_lgammal_r.c:
. Intel 80-bit format implementation of lgammal_r().
* src/e_lgamma.c:
. Expose lgammal as a weak reference to lgamma for platforms
where long double is mapped to double.
* src/e_lgamma_r.c:
. Use integer literal constants instead of real literal constants.
Let compiler(s) do the job of conversion to the appropriate type.
. Expose lgammal_r as a weak reference to lgamma_r for platforms
where long double is mapped to double.
* src/e_lgammaf_r.c:
. Fixed the Cygnus Support conversion of e_lgamma_r.c to float.
This includes the generation of new polynomial and rational
approximations with fewer terms. For each approximation, include
a comment on an estimate of the accuracy over the relevant domain.
. Use integer literal constants instead of real literal constants.
Let compiler(s) do the job of conversion to the appropriate type.
This allows the removal of several explicit casts of double values
to float.
* src/e_lgammal.c:
. Wrapper for lgammal() about lgammal_r().
* src/imprecise.c:
. Remove the lgamma.
* src/math.h:
. Add a prototype for lgammal_r().
* man/lgamma.3:
. Document the new functions.
Reviewed by: bde
Added:
head/lib/msun/ld128/e_lgammal_r.c (contents, props changed)
head/lib/msun/ld80/e_lgammal_r.c (contents, props changed)
head/lib/msun/src/e_lgammal.c (contents, props changed)
Modified:
head/lib/msun/Makefile
head/lib/msun/Symbol.map
head/lib/msun/man/lgamma.3
head/lib/msun/src/e_lgamma.c
head/lib/msun/src/e_lgamma_r.c
head/lib/msun/src/e_lgammaf_r.c
head/lib/msun/src/imprecise.c
head/lib/msun/src/math.h
Modified: head/lib/msun/Makefile
==============================================================================
--- head/lib/msun/Makefile Mon Sep 15 22:32:35 2014 (r271650)
+++ head/lib/msun/Makefile Mon Sep 15 23:21:57 2014 (r271651)
@@ -98,6 +98,7 @@ COMMON_SRCS+= s_copysignl.c s_fabsl.c s_
# If long double != double use these; otherwise, we alias the double versions.
COMMON_SRCS+= e_acoshl.c e_acosl.c e_asinl.c e_atan2l.c e_atanhl.c \
e_coshl.c e_fmodl.c e_hypotl.c \
+ e_lgammal.c e_lgammal_r.c \
e_remainderl.c e_sinhl.c e_sqrtl.c \
invtrig.c k_cosl.c k_sinl.c k_tanl.c \
s_asinhl.c s_atanl.c s_cbrtl.c s_ceill.c s_cosl.c s_cprojl.c \
@@ -188,7 +189,8 @@ MLINKS+=ilogb.3 ilogbf.3 ilogb.3 ilogbl.
ilogb.3 logb.3 ilogb.3 logbf.3 ilogb.3 logbl.3
MLINKS+=j0.3 j1.3 j0.3 jn.3 j0.3 y0.3 j0.3 y1.3 j0.3 y1f.3 j0.3 yn.3
MLINKS+=j0.3 j0f.3 j0.3 j1f.3 j0.3 jnf.3 j0.3 y0f.3 j0.3 ynf.3
-MLINKS+=lgamma.3 gamma.3 lgamma.3 gammaf.3 lgamma.3 lgammaf.3 \
+MLINKS+=lgamma.3 gamma.3 lgamma.3 gammaf.3 \
+ lgamma.3 lgammaf.3 lgamma.3 lgammal.3 \
lgamma.3 tgamma.3 lgamma.3 tgammaf.3
MLINKS+=log.3 log10.3 log.3 log10f.3 log.3 log10l.3 \
log.3 log1p.3 log.3 log1pf.3 log.3 log1pl.3 \
Modified: head/lib/msun/Symbol.map
==============================================================================
--- head/lib/msun/Symbol.map Mon Sep 15 22:32:35 2014 (r271650)
+++ head/lib/msun/Symbol.map Mon Sep 15 23:21:57 2014 (r271651)
@@ -269,6 +269,7 @@ FBSD_1.3 {
erfl;
expl;
expm1l;
+ lgammal;
log10l;
log1pl;
log2l;
@@ -276,7 +277,11 @@ FBSD_1.3 {
sinhl;
tanhl;
/* Implemented as weak aliases for imprecise versions */
- lgammal;
powl;
tgammal;
};
+
+/* First added in 11.0-CURRENT */
+FBSD_1.4 {
+ lgammal_r;
+};
Added: head/lib/msun/ld128/e_lgammal_r.c
==============================================================================
--- /dev/null 00:00:00 1970 (empty, because file is newly added)
+++ head/lib/msun/ld128/e_lgammal_r.c Mon Sep 15 23:21:57 2014
(r271651)
@@ -0,0 +1,329 @@
+/* @(#)e_lgamma_r.c 1.3 95/01/18 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+#include <sys/cdefs.h>
+__FBSDID("$FreeBSD$");
+
+/*
+ * See e_lgamma_r.c for complete comments.
+ *
+ * Converted to long double by Steven G. Kargl.
+ */
+
+#include "fpmath.h"
+#include "math.h"
+#include "math_private.h"
+
+static const volatile double vzero = 0;
+
+static const double
+zero= 0,
+half= 0.5,
+one = 1;
+
+static const long double
+pi = 3.14159265358979323846264338327950288e+00L;
+/*
+ * Domain y in [0x1p-119, 0.28], range ~[-1.4065e-36, 1.4065e-36]:
+ * |(lgamma(2 - y) + y / 2) / y - a(y)| < 2**-119.1
+ */
+static const long double
+a0 = 7.72156649015328606065120900824024296e-02L,
+a1 = 3.22467033424113218236207583323018498e-01L,
+a2 = 6.73523010531980951332460538330282217e-02L,
+a3 = 2.05808084277845478790009252803463129e-02L,
+a4 = 7.38555102867398526627292839296001626e-03L,
+a5 = 2.89051033074152328576829509522483468e-03L,
+a6 = 1.19275391170326097618357349881842913e-03L,
+a7 = 5.09669524743042462515256340206203019e-04L,
+a8 = 2.23154758453578096143609255559576017e-04L,
+a9 = 9.94575127818397632126978731542755129e-05L,
+a10 = 4.49262367375420471287545895027098145e-05L,
+a11 = 2.05072127845117995426519671481628849e-05L,
+a12 = 9.43948816959096748454087141447939513e-06L,
+a13 = 4.37486780697359330303852050718287419e-06L,
+a14 = 2.03920783892362558276037363847651809e-06L,
+a15 = 9.55191070057967287877923073200324649e-07L,
+a16 = 4.48993286185740853170657139487620560e-07L,
+a17 = 2.13107543597620911675316728179563522e-07L,
+a18 = 9.70745379855304499867546549551023473e-08L,
+a19 = 5.61889970390290257926487734695402075e-08L,
+a20 = 6.42739653024130071866684358960960951e-09L,
+a21 = 3.34491062143649291746195612991870119e-08L,
+a22 = -1.57068547394315223934653011440641472e-08L,
+a23 = 1.30812825422415841213733487745200632e-08L;
+/*
+ * Domain x in [tc-0.24, tc+0.28], range ~[-6.3201e-37, 6.3201e-37]:
+ * |(lgamma(x) - tf) - t(x - tc)| < 2**-120.3.
+ */
+static const long double
+tc = 1.46163214496836234126265954232572133e+00L,
+tf = -1.21486290535849608095514557177691584e-01L,
+tt = 1.57061739945077675484237837992951704e-36L,
+t0 = -1.99238329499314692728655623767019240e-36L,
+t1 = -6.08453430711711404116887457663281416e-35L,
+t2 = 4.83836122723810585213722380854828904e-01L,
+t3 = -1.47587722994530702030955093950668275e-01L,
+t4 = 6.46249402389127526561003464202671923e-02L,
+t5 = -3.27885410884813055008502586863748063e-02L,
+t6 = 1.79706751152103942928638276067164935e-02L,
+t7 = -1.03142230366363872751602029672767978e-02L,
+t8 = 6.10053602051788840313573150785080958e-03L,
+t9 = -3.68456960831637325470641021892968954e-03L,
+t10 = 2.25976482322181046611440855340968560e-03L,
+t11 = -1.40225144590445082933490395950664961e-03L,
+t12 = 8.78232634717681264035014878172485575e-04L,
+t13 = -5.54194952796682301220684760591403899e-04L,
+t14 = 3.51912956837848209220421213975000298e-04L,
+t15 = -2.24653443695947456542669289367055542e-04L,
+t16 = 1.44070395420840737695611929680511823e-04L,
+t17 = -9.27609865550394140067059487518862512e-05L,
+t18 = 5.99347334438437081412945428365433073e-05L,
+t19 = -3.88458388854572825603964274134801009e-05L,
+t20 = 2.52476631610328129217896436186551043e-05L,
+t21 = -1.64508584981658692556994212457518536e-05L,
+t22 = 1.07434583475987007495523340296173839e-05L,
+t23 = -7.03070407519397260929482550448878399e-06L,
+t24 = 4.60968590693753579648385629003100469e-06L,
+t25 = -3.02765473778832036018438676945512661e-06L,
+t26 = 1.99238771545503819972741288511303401e-06L,
+t27 = -1.31281299822614084861868817951788579e-06L,
+t28 = 8.60844432267399655055574642052370223e-07L,
+t29 = -5.64535486432397413273248363550536374e-07L,
+t30 = 3.99357783676275660934903139592727737e-07L,
+t31 = -2.95849029193433121795495215869311610e-07L,
+t32 = 1.37790144435073124976696250804940384e-07L;
+/*
+ * Domain y in [-0.1, 0.232], range ~[-1.4046e-37, 1.4181e-37]:
+ * |(lgamma(1 + y) + 0.5 * y) / y - u(y) / v(y)| < 2**-122.8
+ */
+static const long double
+u0 = -7.72156649015328606065120900824024311e-02L,
+u1 = 4.24082772271938167430983113242482656e-01L,
+u2 = 2.96194003481457101058321977413332171e+00L,
+u3 = 6.49503267711258043997790983071543710e+00L,
+u4 = 7.40090051288150177152835698948644483e+00L,
+u5 = 4.94698036296756044610805900340723464e+00L,
+u6 = 2.00194224610796294762469550684947768e+00L,
+u7 = 4.82073087750608895996915051568834949e-01L,
+u8 = 6.46694052280506568192333848437585427e-02L,
+u9 = 4.17685526755100259316625348933108810e-03L,
+u10 = 9.06361003550314327144119307810053410e-05L,
+v1 = 5.15937098592887275994320496999951947e+00L,
+v2 = 1.14068418766251486777604403304717558e+01L,
+v3 = 1.41164839437524744055723871839748489e+01L,
+v4 = 1.07170702656179582805791063277960532e+01L,
+v5 = 5.14448694179047879915042998453632434e+00L,
+v6 = 1.55210088094585540637493826431170289e+00L,
+v7 = 2.82975732849424562719893657416365673e-01L,
+v8 = 2.86424622754753198010525786005443539e-02L,
+v9 = 1.35364253570403771005922441442688978e-03L,
+v10 = 1.91514173702398375346658943749580666e-05L,
+v11 = -3.25364686890242327944584691466034268e-08L;
+/*
+ * Domain x in (2, 3], range ~[-1.3341e-36, 1.3536e-36]:
+ * |(lgamma(y+2) - 0.5 * y) / y - s(y)/r(y)| < 2**-120.1
+ * with y = x - 2.
+ */
+static const long double
+s0 = -7.72156649015328606065120900824024297e-02L,
+s1 = 1.23221687850916448903914170805852253e-01L,
+s2 = 5.43673188699937239808255378293820020e-01L,
+s3 = 6.31998137119005233383666791176301800e-01L,
+s4 = 3.75885340179479850993811501596213763e-01L,
+s5 = 1.31572908743275052623410195011261575e-01L,
+s6 = 2.82528453299138685507186287149699749e-02L,
+s7 = 3.70262021550340817867688714880797019e-03L,
+s8 = 2.83374000312371199625774129290973648e-04L,
+s9 = 1.15091830239148290758883505582343691e-05L,
+s10 = 2.04203474281493971326506384646692446e-07L,
+s11 = 9.79544198078992058548607407635645763e-10L,
+r1 = 2.58037466655605285937112832039537492e+00L,
+r2 = 2.86289413392776399262513849911531180e+00L,
+r3 = 1.78691044735267497452847829579514367e+00L,
+r4 = 6.89400381446725342846854215600008055e-01L,
+r5 = 1.70135865462567955867134197595365343e-01L,
+r6 = 2.68794816183964420375498986152766763e-02L,
+r7 = 2.64617234244861832870088893332006679e-03L,
+r8 = 1.52881761239180800640068128681725702e-04L,
+r9 = 4.63264813762296029824851351257638558e-06L,
+r10 = 5.89461519146957343083848967333671142e-08L,
+r11 = 1.79027678176582527798327441636552968e-10L;
+/*
+ * Domain z in [8, 0x1p70], range ~[-9.8214e-35, 9.8214e-35]:
+ * |lgamma(x) - (x - 0.5) * (log(x) - 1) - w(1/x)| < 2**-113.0
+ */
+static const long double
+w0 = 4.18938533204672741780329736405617738e-01L,
+w1 = 8.33333333333333333333333333332852026e-02L,
+w2 = -2.77777777777777777777777727810123528e-03L,
+w3 = 7.93650793650793650791708939493907380e-04L,
+w4 = -5.95238095238095234390450004444370959e-04L,
+w5 = 8.41750841750837633887817658848845695e-04L,
+w6 = -1.91752691752396849943172337347259743e-03L,
+w7 = 6.41025640880333069429106541459015557e-03L,
+w8 = -2.95506530801732133437990433080327074e-02L,
+w9 = 1.79644237328444101596766586979576927e-01L,
+w10 = -1.39240539108367641920172649259736394e+00L,
+w11 = 1.33987701479007233325288857758641761e+01L,
+w12 = -1.56363596431084279780966590116006255e+02L,
+w13 = 2.14830978044410267201172332952040777e+03L,
+w14 = -3.28636067474227378352761516589092334e+04L,
+w15 = 5.06201257747865138432663574251462485e+05L,
+w16 = -6.79720123352023636706247599728048344e+06L,
+w17 = 6.57556601705472106989497289465949255e+07L,
+w18 = -3.26229058141181783534257632389415580e+08L;
+
+static long double
+sin_pil(long double x)
+{
+ volatile long double vz;
+ long double y,z;
+ uint64_t lx, n;
+ uint16_t hx;
+
+ y = -x;
+
+ vz = y+0x1.p112;
+ z = vz-0x1.p112;
+ if (z == y)
+ return zero;
+
+ vz = y+0x1.p110;
+ EXTRACT_LDBL128_WORDS(hx,lx,n,vz);
+ z = vz-0x1.p110;
+ if (z > y) {
+ z -= 0.25;
+ n--;
+ }
+ n &= 7;
+ y = y - z + n * 0.25L;
+
+ switch (n) {
+ case 0: y = __kernel_sinl(pi*y,zero,0); break;
+ case 1:
+ case 2: y = __kernel_cosl(pi*(0.5-y),zero); break;
+ case 3:
+ case 4: y = __kernel_sinl(pi*(one-y),zero,0); break;
+ case 5:
+ case 6: y = -__kernel_cosl(pi*(y-1.5),zero); break;
+ default: y = __kernel_sinl(pi*(y-2.0),zero,0); break;
+ }
+ return -y;
+}
+
+
+long double
+lgammal_r(long double x, int *signgamp)
+{
+ long double nadj,p,p1,p2,p3,q,r,t,w,y,z;
+ uint64_t llx,lx;
+ int i;
+ uint16_t hx;
+
+ EXTRACT_LDBL128_WORDS(hx, lx, llx, x);
+
+ if((hx & 0x7fff) == 0x7fff) { /* erfl(nan)=nan */
+ i = (hx>>15)<<1;
+ return (1-i)+one/x; /* erfl(+-inf)=+-1 */
+ }
+
+ /* purge off +-inf, NaN, +-0, tiny and negative arguments */
+ *signgamp = 1;
+ if((hx & 0x7fff) == 0x7fff) /* x is +-Inf or NaN */
+ return x*x;
+ if((hx==0||hx==0x8000)&&lx==0) return one/vzero;
+
+ /* purge off tiny and negative arguments */
+ if(fabsl(x)<0x1p-119L) {
+ if(hx&0x8000) {
+ *signgamp = -1;
+ return -logl(-x);
+ } else return -logl(x);
+ }
+ if(hx&0x8000) {
+ if(fabsl(x)>=0x1p112)
+ return one/vzero;
+ t = sin_pil(x);
+ if(t==zero) return one/vzero;
+ nadj = logl(pi/fabsl(t*x));
+ if(t<zero) *signgamp = -1;
+ x = -x;
+ }
+
+ if(x == 1 || x ==2) r = 0;
+ else if(x<2) {
+ if(x<=0.8999996185302734) {
+ r = -logl(x);
+ if(x>=0.7315998077392578) {y = 1-x; i= 0;}
+ else if(x>=0.2316399812698364) {y= x-(tc-1); i=1;}
+ else {y = x; i=2;}
+ } else {
+ r = 0;
+ if(x>=1.7316312789916992) {y=2-x;i=0;}
+ else if(x>=1.2316322326660156) {y=x-tc;i=1;}
+ else {y=x-1;i=2;}
+ }
+ switch(i) {
+ case 0:
+ z = y*y;
+ p1 = a0+z*(a2+z*(a4+z*(a6+z*(a8+z*(a10+z*(a12+z*(a14+z*(a16+
+ z*(a18+z*(a20+z*a22))))))))));
+ p2 = z*(a1+z*(a3+z*(a5+z*(a7+z*(a9+z*(a11+z*(a13+z*(a15+
+ z*(a17+z*(a19+z*(a21+z*a23)))))))))));
+ p = y*p1+p2;
+ r += (p-y/2); break;
+ case 1:
+ p = t0+y*t1+tt+y*y*(t2+y*(t3+y*(t4+y*(t5+y*(t6+y*(t7+y*(t8+
+ y*(t9+y*(t10+y*(t11+y*(t12+y*(t13+y*(t14+y*(t15+y*(t16+
+ y*(t17+y*(t18+y*(t19+y*(t20+y*(t21+y*(t22+y*(t23+
+ y*(t24+y*(t25+y*(t26+y*(t27+y*(t28+y*(t29+y*(t30+
+ y*(t31+y*t32))))))))))))))))))))))))))))));
+ r += (tf + p); break;
+ case 2:
+ p1 = y*(u0+y*(u1+y*(u2+y*(u3+y*(u4+y*(u5+y*(u6+y*(u7+
+ y*(u8+y*(u9+y*u10))))))))));
+ p2 = one+y*(v1+y*(v2+y*(v3+y*(v4+y*(v5+y*(v6+y*(v7+
+ y*(v8+y*(v9+y*(v10+y*v11))))))))));
+ r += (-y/2 + p1/p2);
+ }
+ }
+ else if(x<8) {
+ i = x;
+ y = x-i;
+ p = y*(s0+y*(s1+y*(s2+y*(s3+y*(s4+y*(s5+y*(s6+y*(s7+y*(s8+
+ y*(s9+y*(s10+y*s11)))))))))));
+ q = one+y*(r1+y*(r2+y*(r3+y*(r4+y*(r5+y*(r6+y*(r7+y*(r8+
+ y*(r9+y*(r10+y*r11))))))))));
+ r = y/2+p/q;
+ z = 1; /* lgamma(1+s) = log(s) + lgamma(s) */
+ switch(i) {
+ case 7: z *= (y+6); /* FALLTHRU */
+ case 6: z *= (y+5); /* FALLTHRU */
+ case 5: z *= (y+4); /* FALLTHRU */
+ case 4: z *= (y+3); /* FALLTHRU */
+ case 3: z *= (y+2); /* FALLTHRU */
+ r += logl(z); break;
+ }
+ } else if (x < 0x1p119L) {
+ t = logl(x);
+ z = one/x;
+ y = z*z;
+ w = w0+z*(w1+y*(w2+y*(w3+y*(w4+y*(w5+y*(w6+y*(w7+y*(w8+
+ y*(w9+y*(w10+y*(w11+y*(w12+y*(w13+y*(w14+y*(w15+y*(w16+
+ y*(w17+y*w18)))))))))))))))));
+ r = (x-half)*(t-one)+w;
+ } else
+ r = x*(logl(x)-1);
+ if(hx&0x8000) r = nadj - r;
+ return r;
+}
Added: head/lib/msun/ld80/e_lgammal_r.c
==============================================================================
--- /dev/null 00:00:00 1970 (empty, because file is newly added)
+++ head/lib/msun/ld80/e_lgammal_r.c Mon Sep 15 23:21:57 2014
(r271651)
@@ -0,0 +1,345 @@
+/* @(#)e_lgamma_r.c 1.3 95/01/18 */
+/*
+ * ====================================================
+ * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+#include <sys/cdefs.h>
+__FBSDID("$FreeBSD$");
+
+/*
+ * See s_lgamma_r.c for complete comments.
+ *
+ * Converted to long double by Steven G. Kargl.
+ */
+
+#include <float.h>
+#ifdef __i386__
+#include <ieeefp.h>
+#endif
+
+#include "fpmath.h"
+#include "math.h"
+#include "math_private.h"
+
+static const volatile double vzero = 0;
+
+static const double
+zero= 0,
+half= 0.5,
+one = 1;
+
+static const union IEEEl2bits
+piu = LD80C(0xc90fdaa22168c235, 1, 3.14159265358979323851e+00L);
+#define pi (piu.e)
+/*
+ * Domain y in [0x1p-70, 0.27], range ~[-4.5264e-22, 4.5264e-22]:
+ * |(lgamma(2 - y) + y / 2) / y - a(y)| < 2**-70.9
+ */
+static const union IEEEl2bits
+a0u = LD80C(0x9e233f1bed863d26, -4, 7.72156649015328606028e-02L),
+a1u = LD80C(0xa51a6625307d3249, -2, 3.22467033424113218889e-01L),
+a2u = LD80C(0x89f000d2abafda8c, -4, 6.73523010531979398946e-02L),
+a3u = LD80C(0xa8991563eca75f26, -6, 2.05808084277991211934e-02L),
+a4u = LD80C(0xf2027e10634ce6b6, -8, 7.38555102796070454026e-03L),
+a5u = LD80C(0xbd6eb76dd22187f4, -9, 2.89051035162703932972e-03L),
+a6u = LD80C(0x9c562ab05e0458ed, -10, 1.19275351624639999297e-03L),
+a7u = LD80C(0x859baed93ee48e46, -11, 5.09674593842117925320e-04L),
+a8u = LD80C(0xe9f28a4432949af2, -13, 2.23109648015769155122e-04L),
+a9u = LD80C(0xd12ad0d9b93c6bb0, -14, 9.97387167479808509830e-05L),
+a10u= LD80C(0xb7522643c78a219b, -15, 4.37071076331030136818e-05L),
+a11u= LD80C(0xca024dcdece2cb79, -16, 2.40813493372040143061e-05L),
+a12u= LD80C(0xbb90fb6968ebdbf9, -19, 2.79495621083634031729e-06L),
+a13u= LD80C(0xba1c9ffeeae07b37, -17, 1.10931287015513924136e-05L);
+#define a0 (a0u.e)
+#define a1 (a1u.e)
+#define a2 (a2u.e)
+#define a3 (a3u.e)
+#define a4 (a4u.e)
+#define a5 (a5u.e)
+#define a6 (a6u.e)
+#define a7 (a7u.e)
+#define a8 (a8u.e)
+#define a9 (a9u.e)
+#define a10 (a10u.e)
+#define a11 (a11u.e)
+#define a12 (a12u.e)
+#define a13 (a13u.e)
+/*
+ * Domain x in [tc-0.24, tc+0.28], range ~[-6.1205e-22, 6.1205e-22]:
+ * |(lgamma(x) - tf) - t(x - tc)| < 2**-70.5
+ */
+static const union IEEEl2bits
+tcu = LD80C(0xbb16c31ab5f1fb71, 0, 1.46163214496836234128e+00L),
+tfu = LD80C(0xf8cdcde61c520e0f, -4, -1.21486290535849608093e-01L),
+ttu = LD80C(0xd46ee54b27d4de99, -69, -2.81152980996018785880e-21L),
+t0u = LD80C(0x80b9406556a62a6b, -68, 3.40728634996055147231e-21L),
+t1u = LD80C(0xc7e9c6f6df3f8c39, -67, -1.05833162742737073665e-20L),
+t2u = LD80C(0xf7b95e4771c55d51, -2, 4.83836122723810583532e-01L),
+t3u = LD80C(0x97213c6e35e119ff, -3, -1.47587722994530691476e-01L),
+t4u = LD80C(0x845a14a6a81dc94b, -4, 6.46249402389135358063e-02L),
+t5u = LD80C(0x864d46fa89997796, -5, -3.27885410884846056084e-02L),
+t6u = LD80C(0x93373cbd00297438, -6, 1.79706751150707171293e-02L),
+t7u = LD80C(0xa8fcfca7eddc8d1d, -7, -1.03142230361450732547e-02L),
+t8u = LD80C(0xc7e7015ff4bc45af, -8, 6.10053603296546099193e-03L),
+t9u = LD80C(0xf178d2247adc5093, -9, -3.68456964904901200152e-03L),
+t10u = LD80C(0x94188d58f12e5e9f, -9, 2.25976420273774583089e-03L),
+t11u = LD80C(0xb7cbaef14e1406f1, -10, -1.40224943666225639823e-03L),
+t12u = LD80C(0xe63a671e6704ea4d, -11, 8.78250640744776944887e-04L),
+t13u = LD80C(0x914b6c9cae61783e, -11, -5.54255012657716808811e-04L),
+t14u = LD80C(0xb858f5bdb79276fe, -12, 3.51614951536825927370e-04L),
+t15u = LD80C(0xea73e744c34b9591, -13, -2.23591563824520112236e-04L),
+t16u = LD80C(0x99aeabb0d67ba835, -13, 1.46562869351659194136e-04L),
+t17u = LD80C(0xd7c6938325db2024, -14, -1.02889866046435680588e-04L),
+t18u = LD80C(0xe24cb1e3b0474775, -15, 5.39540265505221957652e-05L);
+#define tc (tcu.e)
+#define tf (tfu.e)
+#define tt (ttu.e)
+#define t0 (t0u.e)
+#define t1 (t1u.e)
+#define t2 (t2u.e)
+#define t3 (t3u.e)
+#define t4 (t4u.e)
+#define t5 (t5u.e)
+#define t6 (t6u.e)
+#define t7 (t7u.e)
+#define t8 (t8u.e)
+#define t9 (t9u.e)
+#define t10 (t10u.e)
+#define t11 (t11u.e)
+#define t12 (t12u.e)
+#define t13 (t13u.e)
+#define t14 (t14u.e)
+#define t15 (t15u.e)
+#define t16 (t16u.e)
+#define t17 (t17u.e)
+#define t18 (t18u.e)
+/*
+ * Domain y in [-0.1, 0.232], range ~[-8.1938e-22, 8.3815e-22]:
+ * |(lgamma(1 + y) + 0.5 * y) / y - u(y) / v(y)| < 2**-71.2
+ */
+static const union IEEEl2bits
+u0u = LD80C(0x9e233f1bed863d27, -4, -7.72156649015328606095e-02L),
+u1u = LD80C(0x98280ee45e4ddd3d, -1, 5.94361239198682739769e-01L),
+u2u = LD80C(0xe330c8ead4130733, 0, 1.77492629495841234275e+00L),
+u3u = LD80C(0xd4a213f1a002ec52, 0, 1.66119622514818078064e+00L),
+u4u = LD80C(0xa5a9ca6f5bc62163, -1, 6.47122051417476492989e-01L),
+u5u = LD80C(0xc980e49cd5b019e6, -4, 9.83903751718671509455e-02L),
+u6u = LD80C(0xff636a8bdce7025b, -9, 3.89691687802305743450e-03L),
+v1u = LD80C(0xbd109c533a19fbf5, 1, 2.95413883330948556544e+00L),
+v2u = LD80C(0xd295cbf96f31f099, 1, 3.29039286955665403176e+00L),
+v3u = LD80C(0xdab8bcfee40496cb, 0, 1.70876276441416471410e+00L),
+v4u = LD80C(0xd2f2dc3638567e9f, -2, 4.12009126299534668571e-01L),
+v5u = LD80C(0xa07d9b0851070f41, -5, 3.91822868305682491442e-02L),
+v6u = LD80C(0xe3cd8318f7adb2c4, -11, 8.68998648222144351114e-04L);
+#define u0 (u0u.e)
+#define u1 (u1u.e)
+#define u2 (u2u.e)
+#define u3 (u3u.e)
+#define u4 (u4u.e)
+#define u5 (u5u.e)
+#define u6 (u6u.e)
+#define v1 (v1u.e)
+#define v2 (v2u.e)
+#define v3 (v3u.e)
+#define v4 (v4u.e)
+#define v5 (v5u.e)
+#define v6 (v6u.e)
+/*
+ * Domain x in (2, 3], range ~[-3.3648e-22, 3.4416e-22]:
+ * |(lgamma(y+2) - 0.5 * y) / y - s(y)/r(y)| < 2**-72.3
+ * with y = x - 2.
+ */
+static const union IEEEl2bits
+s0u = LD80C(0x9e233f1bed863d27, -4, -7.72156649015328606095e-02L),
+s1u = LD80C(0xd3ff0dcc7fa91f94, -3, 2.07027640921219389860e-01L),
+s2u = LD80C(0xb2bb62782478ef31, -2, 3.49085881391362090549e-01L),
+s3u = LD80C(0xb49f7438c4611a74, -3, 1.76389518704213357954e-01L),
+s4u = LD80C(0x9a957008fa27ecf9, -5, 3.77401710862930008071e-02L),
+s5u = LD80C(0xda9b389a6ca7a7ac, -9, 3.33566791452943399399e-03L),
+s6u = LD80C(0xbc7a2263faf59c14, -14, 8.98728786745638844395e-05L),
+r1u = LD80C(0xbf5cff5b11477d4d, 0, 1.49502555796294337722e+00L),
+r2u = LD80C(0xd9aec89de08e3da6, -1, 8.50323236984473285866e-01L),
+r3u = LD80C(0xeab7ae5057c443f9, -3, 2.29216312078225806131e-01L),
+r4u = LD80C(0xf29707d9bd2b1e37, -6, 2.96130326586640089145e-02L),
+r5u = LD80C(0xd376c2f09736c5a3, -10, 1.61334161411590662495e-03L),
+r6u = LD80C(0xc985983d0cd34e3d, -16, 2.40232770710953450636e-05L),
+r7u = LD80C(0xe5c7a4f7fc2ef13d, -25, -5.34997929289167573510e-08L);
+#define s0 (s0u.e)
+#define s1 (s1u.e)
+#define s2 (s2u.e)
+#define s3 (s3u.e)
+#define s4 (s4u.e)
+#define s5 (s5u.e)
+#define s6 (s6u.e)
+#define r1 (r1u.e)
+#define r2 (r2u.e)
+#define r3 (r3u.e)
+#define r4 (r4u.e)
+#define r5 (r5u.e)
+#define r6 (r6u.e)
+#define r7 (r7u.e)
+/*
+ * Domain z in [8, 0x1p70], range ~[-3.0235e-22, 3.0563e-22]:
+ * |lgamma(x) - (x - 0.5) * (log(x) - 1) - w(1/x)| < 2**-71.7
+ */
+static const union IEEEl2bits
+w0u = LD80C(0xd67f1c864beb4a69, -2, 4.18938533204672741776e-01L),
+w1u = LD80C(0xaaaaaaaaaaaaaaa1, -4, 8.33333333333333332678e-02L),
+w2u = LD80C(0xb60b60b60b5491c9, -9, -2.77777777777760927870e-03L),
+w3u = LD80C(0xd00d00cf58aede4c, -11, 7.93650793490637233668e-04L),
+w4u = LD80C(0x9c09bf626783d4a5, -11, -5.95238023926039051268e-04L),
+w5u = LD80C(0xdca7cadc5baa517b, -11, 8.41733700408000822962e-04L),
+w6u = LD80C(0xfb060e361e1ffd07, -10, -1.91515849570245136604e-03L),
+w7u = LD80C(0xcbd5101bb58d1f2b, -8, 6.22046743903262649294e-03L),
+w8u = LD80C(0xad27a668d32c821b, -6, -2.11370706734662081843e-02L);
+#define w0 (w0u.e)
+#define w1 (w1u.e)
+#define w2 (w2u.e)
+#define w3 (w3u.e)
+#define w4 (w4u.e)
+#define w5 (w5u.e)
+#define w6 (w6u.e)
+#define w7 (w7u.e)
+#define w8 (w8u.e)
+
+static long double
+sin_pil(long double x)
+{
+ volatile long double vz;
+ long double y,z;
+ uint64_t n;
+ uint16_t hx;
+
+ y = -x;
+
+ vz = y+0x1p63L;
+ z = vz-0x1p63L;
+ if (z == y)
+ return zero;
+
+ vz = y+0x1p61;
+ EXTRACT_LDBL80_WORDS(hx,n,vz);
+ z = vz-0x1p61;
+ if (z > y) {
+ z -= 0.25; /* adjust to round down */
+ n--;
+ }
+ n &= 7; /* octant of y mod 2 */
+ y = y - z + n * 0.25; /* y mod 2 */
+
+ switch (n) {
+ case 0: y = __kernel_sinl(pi*y,zero,0); break;
+ case 1:
+ case 2: y = __kernel_cosl(pi*(0.5-y),zero); break;
+ case 3:
+ case 4: y = __kernel_sinl(pi*(one-y),zero,0); break;
+ case 5:
+ case 6: y = -__kernel_cosl(pi*(y-1.5),zero); break;
+ default: y = __kernel_sinl(pi*(y-2.0),zero,0); break;
+ }
+ return -y;
+}
+
+long double
+lgammal_r(long double x, int *signgamp)
+{
+ long double nadj,p,p1,p2,p3,q,r,t,w,y,z;
+ uint64_t lx;
+ int i;
+ uint16_t hx;
+
+ EXTRACT_LDBL80_WORDS(hx,lx,x);
+
+ /* purge off +-inf, NaN, +-0 */
+ *signgamp = 1;
+ if((hx & 0x7fff) == 0x7fff) /* x is +-Inf or NaN */
+ return x*x;
+ if((hx==0||hx==0x8000)&&lx==0) return one/vzero;
+
+ ENTERI();
+
+ /* purge off tiny and negative arguments */
+ if(fabsl(x)<0x1p-70L) { /* |x|<2**-70, return -log(|x|) */
+ if(hx&0x8000) {
+ *signgamp = -1;
+ RETURNI(-logl(-x));
+ } else RETURNI(-logl(x));
+ }
+ if(hx&0x8000) {
+ if(fabsl(x)>=0x1p63) /* |x|>=2**(p-1), must be -integer */
+ RETURNI(one/vzero);
+ t = sin_pil(x);
+ if(t==zero) RETURNI(one/vzero); /* -integer */
+ nadj = logl(pi/fabsl(t*x));
+ if(t<zero) *signgamp = -1;
+ x = -x;
+ }
+
+ /* purge off 1 and 2 */
+ if(x == 1 || x == 2) r = 0;
+ /* for x < 2.0 */
+ else if(x<2) {
+ if(x<=0.8999996185302734) { /* lgamma(x) = lgamma(x+1)-log(x) */
+ r = - logl(x);
+ if(x>=0.7315998077392578) {y = 1-x; i= 0;}
+ else if(x>=0.2316399812698364) {y= x-(tc-1); i=1;}
+ else {y = x; i=2;}
+ } else {
+ r = 0;
+ if(x>=1.7316312789916992) {y=2-x;i=0;}
+ else if(x>=1.2316322326660156) {y=x-tc;i=1;}
+ else {y=x-1;i=2;}
+ }
+ switch(i) {
+ case 0:
+ z = y*y;
+ p1 = a0+z*(a2+z*(a4+z*(a6+z*(a8+z*(a10+z*a12)))));
+ p2 = z*(a1+z*(a3+z*(a5+z*(a7+z*(a9+z*(a11+z*a13))))));
+ p = y*p1+p2;
+ r += (p-y/2); break;
+ case 1:
+ p = t0+y*t1+tt+y*y*(t2+y*(t3+y*(t4+y*(t5+y*(t6+y*(t7+y*(t8+
+ y*(t9+y*(t10+y*(t11+y*(t12+y*(t13+y*(t14+y*(t15+y*(t16+
+ y*(t17+y*t18))))))))))))))));
+ r += (tf + p); break;
+ case 2:
+ p1 = y*(u0+y*(u1+y*(u2+y*(u3+y*(u4+y*(u5+y*u6))))));
+ p2 = 1+y*(v1+y*(v2+y*(v3+y*(v4+y*(v5+y*v6)))));
+ r += (-y/2 + p1/p2);
+ }
+ }
+ else if(x<8) {
+ i = x;
+ y = x-i;
+ p = y*(s0+y*(s1+y*(s2+y*(s3+y*(s4+y*(s5+y*s6))))));
+ q = 1+y*(r1+y*(r2+y*(r3+y*(r4+y*(r5+y*(r6+y*r7))))));
+ r = y/2+p/q;
+ z = 1; /* lgamma(1+s) = log(s) + lgamma(s) */
+ switch(i) {
+ case 7: z *= (y+6); /* FALLTHRU */
+ case 6: z *= (y+5); /* FALLTHRU */
+ case 5: z *= (y+4); /* FALLTHRU */
+ case 4: z *= (y+3); /* FALLTHRU */
+ case 3: z *= (y+2); /* FALLTHRU */
+ r += logl(z); break;
+ }
+ /* 8.0 <= x < 2**70 */
+ } else if (x < 0x1p70L) {
+ t = logl(x);
+ z = one/x;
+ y = z*z;
+ w = w0+z*(w1+y*(w2+y*(w3+y*(w4+y*(w5+y*(w6+y*(w7+y*w8)))))));
+ r = (x-half)*(t-one)+w;
+ } else
+ /* 2**70 <= x <= inf */
+ r = x*(logl(x)-1);
+ if(hx&0x8000) r = nadj - r;
+ RETURNI(r);
+}
Modified: head/lib/msun/man/lgamma.3
==============================================================================
--- head/lib/msun/man/lgamma.3 Mon Sep 15 22:32:35 2014 (r271650)
+++ head/lib/msun/man/lgamma.3 Mon Sep 15 23:21:57 2014 (r271651)
@@ -28,7 +28,7 @@
.\" from: @(#)lgamma.3 6.6 (Berkeley) 12/3/92
.\" $FreeBSD$
.\"
-.Dd January 14, 2005
+.Dd September 12, 2014
.Dt LGAMMA 3
.Os
.Sh NAME
@@ -36,6 +36,8 @@
.Nm lgamma_r ,
.Nm lgammaf ,
.Nm lgammaf_r ,
+.Nm lgammal ,
+.Nm lgammal_r ,
.Nm gamma ,
.Nm gamma_r ,
.Nm gammaf ,
@@ -58,6 +60,10 @@
.Fn lgammaf "float x"
.Ft float
.Fn lgammaf_r "float x" "int *signgamp"
+.Ft "long double"
+.Fn lgammal "long double x"
+.Ft "long double"
+.Fn lgammal_r "long double x" "int *signgamp"
.Ft double
.Fn gamma "double x"
.Ft double
@@ -66,14 +72,15 @@
.Fn gammaf "float x"
.Ft float
.Fn gammaf_r "float x" "int *signgamp"
-.Ft double
+.Ft "long double"
.Fn tgamma "double x"
.Ft float
.Fn tgammaf "float x"
.Sh DESCRIPTION
-.Fn lgamma x
+.Fn lgamma x ,
+.Fn lgammaf x ,
and
-.Fn lgammaf x
+.Fn lgammal x
.if t \{\
return ln\||\(*G(x)| where
.Bd -unfilled -offset indent
@@ -87,13 +94,15 @@ The external integer
.Fa signgam
returns the sign of \(*G(x).
.Pp
-.Fn lgamma_r x signgamp
+.Fn lgamma_r x signgamp ,
+.Fn lgammaf_r x signgamp ,
and
-.Fn lgammaf_r x signgamp
+.Fn lgammal_r x signgamp
provide the same functionality as
-.Fn lgamma x
+.Fn lgamma x ,
+.Fn lgammaf x ,
and
-.Fn lgammaf x
+.Fn lgammal x ,
but the caller must provide an integer to store the sign of \(*G(x).
.Pp
The
@@ -115,6 +124,7 @@ are deprecated aliases for
and
.Fn lgammaf_r ,
respectively.
+
.Sh IDIOSYNCRASIES
Do not use the expression
.Dq Li signgam\(**exp(lgamma(x))
@@ -139,14 +149,18 @@ Exponentiation of
will lose up to 10 significant bits.
.Sh RETURN VALUES
.Fn gamma ,
-.Fn gamma_r ,
.Fn gammaf ,
+.Fn gammal ,
+.Fn gamma_r ,
.Fn gammaf_r ,
+.Fn gammal_r ,
.Fn lgamma ,
-.Fn lgamma_r ,
.Fn lgammaf ,
+.Fn lgammal ,
+.Fn lgamma_r ,
+.Fn lgammaf_r ,
and
-.Fn lgammaf_r
+.Fn lgammal_r
return appropriate values unless an argument is out of range.
Overflow will occur for sufficiently large positive values, and
non-positive integers.
@@ -159,6 +173,7 @@ will underflow.
The
.Fn lgamma ,
.Fn lgammaf ,
+.Fn lgammal ,
.Fn tgamma ,
and
.Fn tgammaf
Modified: head/lib/msun/src/e_lgamma.c
==============================================================================
--- head/lib/msun/src/e_lgamma.c Mon Sep 15 22:32:35 2014
(r271650)
+++ head/lib/msun/src/e_lgamma.c Mon Sep 15 23:21:57 2014
(r271651)
@@ -21,6 +21,8 @@ __FBSDID("$FreeBSD$");
* Method: call __ieee754_lgamma_r
*/
+#include <float.h>
+
#include "math.h"
#include "math_private.h"
@@ -31,3 +33,7 @@ __ieee754_lgamma(double x)
{
return __ieee754_lgamma_r(x,&signgam);
}
+
+#if (LDBL_MANT_DIG == 53)
+__weak_reference(lgamma, lgammal);
+#endif
Modified: head/lib/msun/src/e_lgamma_r.c
==============================================================================
--- head/lib/msun/src/e_lgamma_r.c Mon Sep 15 22:32:35 2014
(r271650)
+++ head/lib/msun/src/e_lgamma_r.c Mon Sep 15 23:21:57 2014
(r271651)
@@ -83,6 +83,8 @@ __FBSDID("$FreeBSD$");
*
*/
+#include <float.h>
+
#include "math.h"
#include "math_private.h"
@@ -250,7 +252,7 @@ __ieee754_lgamma_r(double x, int *signga
p1 = a0+z*(a2+z*(a4+z*(a6+z*(a8+z*a10))));
p2 = z*(a1+z*(a3+z*(a5+z*(a7+z*(a9+z*a11)))));
p = y*p1+p2;
- r += (p-0.5*y); break;
+ r += (p-y/2); break;
case 1:
z = y*y;
w = z*y;
@@ -273,11 +275,11 @@ __ieee754_lgamma_r(double x, int *signga
r = half*y+p/q;
z = one; /* lgamma(1+s) = log(s) + lgamma(s) */
switch(i) {
- case 7: z *= (y+6.0); /* FALLTHRU */
- case 6: z *= (y+5.0); /* FALLTHRU */
- case 5: z *= (y+4.0); /* FALLTHRU */
- case 4: z *= (y+3.0); /* FALLTHRU */
- case 3: z *= (y+2.0); /* FALLTHRU */
+ case 7: z *= (y+6); /* FALLTHRU */
+ case 6: z *= (y+5); /* FALLTHRU */
+ case 5: z *= (y+4); /* FALLTHRU */
+ case 4: z *= (y+3); /* FALLTHRU */
+ case 3: z *= (y+2); /* FALLTHRU */
r += __ieee754_log(z); break;
}
/* 8.0 <= x < 2**58 */
@@ -293,3 +295,8 @@ __ieee754_lgamma_r(double x, int *signga
if(hx<0) r = nadj - r;
return r;
}
+
+#if (LDBL_MANT_DIG == 53)
+__weak_reference(lgamma_r, lgammal_r);
+#endif
+
Modified: head/lib/msun/src/e_lgammaf_r.c
==============================================================================
--- head/lib/msun/src/e_lgammaf_r.c Mon Sep 15 22:32:35 2014
(r271650)
+++ head/lib/msun/src/e_lgammaf_r.c Mon Sep 15 23:21:57 2014
(r271651)
@@ -1,5 +1,6 @@
/* e_lgammaf_r.c -- float version of e_lgamma_r.c.
* Conversion to float by Ian Lance Taylor, Cygnus Support, i...@cygnus.com.
+ * Conversion to float fixed By Steven G. Kargl.
*/
/*
@@ -22,72 +23,63 @@ __FBSDID("$FreeBSD$");
static const volatile float vzero = 0;
static const float
-zero= 0.0000000000e+00,
-half= 5.0000000000e-01, /* 0x3f000000 */
-one = 1.0000000000e+00, /* 0x3f800000 */
+zero= 0,
+half= 0.5,
+one = 1,
pi = 3.1415927410e+00, /* 0x40490fdb */
-a0 = 7.7215664089e-02, /* 0x3d9e233f */
-a1 = 3.2246702909e-01, /* 0x3ea51a66 */
-a2 = 6.7352302372e-02, /* 0x3d89f001 */
-a3 = 2.0580807701e-02, /* 0x3ca89915 */
-a4 = 7.3855509982e-03, /* 0x3bf2027e */
-a5 = 2.8905137442e-03, /* 0x3b3d6ec6 */
-a6 = 1.1927076848e-03, /* 0x3a9c54a1 */
-a7 = 5.1006977446e-04, /* 0x3a05b634 */
-a8 = 2.2086278477e-04, /* 0x39679767 */
-a9 = 1.0801156895e-04, /* 0x38e28445 */
-a10 = 2.5214456400e-05, /* 0x37d383a2 */
-a11 = 4.4864096708e-05, /* 0x383c2c75 */
-tc = 1.4616321325e+00, /* 0x3fbb16c3 */
-tf = -1.2148628384e-01, /* 0xbdf8cdcd */
-/* tt = -(tail of tf) */
-tt = 6.6971006518e-09, /* 0x31e61c52 */
-t0 = 4.8383611441e-01, /* 0x3ef7b95e */
-t1 = -1.4758771658e-01, /* 0xbe17213c */
-t2 = 6.4624942839e-02, /* 0x3d845a15 */
-t3 = -3.2788541168e-02, /* 0xbd064d47 */
-t4 = 1.7970675603e-02, /* 0x3c93373d */
-t5 = -1.0314224288e-02, /* 0xbc28fcfe */
-t6 = 6.1005386524e-03, /* 0x3bc7e707 */
-t7 = -3.6845202558e-03, /* 0xbb7177fe */
-t8 = 2.2596477065e-03, /* 0x3b141699 */
-t9 = -1.4034647029e-03, /* 0xbab7f476 */
-t10 = 8.8108185446e-04, /* 0x3a66f867 */
-t11 = -5.3859531181e-04, /* 0xba0d3085 */
-t12 = 3.1563205994e-04, /* 0x39a57b6b */
-t13 = -3.1275415677e-04, /* 0xb9a3f927 */
-t14 = 3.3552918467e-04, /* 0x39afe9f7 */
-u0 = -7.7215664089e-02, /* 0xbd9e233f */
-u1 = 6.3282704353e-01, /* 0x3f2200f4 */
-u2 = 1.4549225569e+00, /* 0x3fba3ae7 */
-u3 = 9.7771751881e-01, /* 0x3f7a4bb2 */
-u4 = 2.2896373272e-01, /* 0x3e6a7578 */
-u5 = 1.3381091878e-02, /* 0x3c5b3c5e */
-v1 = 2.4559779167e+00, /* 0x401d2ebe */
-v2 = 2.1284897327e+00, /* 0x4008392d */
-v3 = 7.6928514242e-01, /* 0x3f44efdf */
-v4 = 1.0422264785e-01, /* 0x3dd572af */
-v5 = 3.2170924824e-03, /* 0x3b52d5db */
-s0 = -7.7215664089e-02, /* 0xbd9e233f */
-s1 = 2.1498242021e-01, /* 0x3e5c245a */
-s2 = 3.2577878237e-01, /* 0x3ea6cc7a */
-s3 = 1.4635047317e-01, /* 0x3e15dce6 */
-s4 = 2.6642270386e-02, /* 0x3cda40e4 */
-s5 = 1.8402845599e-03, /* 0x3af135b4 */
-s6 = 3.1947532989e-05, /* 0x3805ff67 */
-r1 = 1.3920053244e+00, /* 0x3fb22d3b */
-r2 = 7.2193557024e-01, /* 0x3f38d0c5 */
-r3 = 1.7193385959e-01, /* 0x3e300f6e */
-r4 = 1.8645919859e-02, /* 0x3c98bf54 */
-r5 = 7.7794247773e-04, /* 0x3a4beed6 */
-r6 = 7.3266842264e-06, /* 0x36f5d7bd */
-w0 = 4.1893854737e-01, /* 0x3ed67f1d */
-w1 = 8.3333335817e-02, /* 0x3daaaaab */
-w2 = -2.7777778450e-03, /* 0xbb360b61 */
-w3 = 7.9365057172e-04, /* 0x3a500cfd */
-w4 = -5.9518753551e-04, /* 0xba1c065c */
-w5 = 8.3633989561e-04, /* 0x3a5b3dd2 */
-w6 = -1.6309292987e-03; /* 0xbad5c4e8 */
+/*
+ * Domain y in [0x1p-27, 0.27], range ~[-3.4599e-10, 3.4590e-10]:
+ * |(lgamma(2 - y) + 0.5 * y) / y - a(y)| < 2**-31.4
+ */
+a0 = 7.72156641e-02, /* 0x3d9e233f */
+a1 = 3.22467119e-01, /* 0x3ea51a69 */
+a2 = 6.73484802e-02, /* 0x3d89ee00 */
+a3 = 2.06395667e-02, /* 0x3ca9144f */
+a4 = 6.98275631e-03, /* 0x3be4cf9b */
+a5 = 4.11768444e-03, /* 0x3b86eda4 */
+/*
+ * Domain x in [tc-0.24, tc+0.28], range ~[-5.6577e-10, 5.5677e-10]:
+ * |(lgamma(x) - tf) - t(x - tc)| < 2**-30.8.
+ */
+tc = 1.46163213e+00, /* 0x3fbb16c3 */
+tf = -1.21486291e-01, /* 0xbdf8cdce */
+t0 = -2.94064460e-11, /* 0xae0154b7 */
+t1 = -2.35939837e-08, /* 0xb2caabb8 */
+t2 = 4.83836412e-01, /* 0x3ef7b968 */
+t3 = -1.47586212e-01, /* 0xbe1720d7 */
+t4 = 6.46013096e-02, /* 0x3d844db1 */
+t5 = -3.28450352e-02, /* 0xbd068884 */
*** DIFF OUTPUT TRUNCATED AT 1000 LINES ***
_______________________________________________
svn-src-head@freebsd.org mailing list
http://lists.freebsd.org/mailman/listinfo/svn-src-head
To unsubscribe, send any mail to "svn-src-head-unsubscr...@freebsd.org"
