summaryrefslogtreecommitdiff
path: root/apps
diff options
context:
space:
mode:
authorFranklin Wei <frankhwei536@gmail.com>2016-11-27 20:55:46 -0500
committerFranklin Wei <frankhwei536@gmail.com>2016-11-27 20:55:46 -0500
commitbad4eb382b571e3b883208980719bcc065218712 (patch)
tree20fde33eab8628195c95c7bf6c573fe653260f25 /apps
parent3993f11a59fd912dba3a789f2c17fd6073b02688 (diff)
parentd71a2abe90a751fc1b8321fc53c38098c17b6bc8 (diff)
downloadrockbox-working.zip
rockbox-working.tar.gz
rockbox-working.tar.bz2
rockbox-working.tar.xz
Merge branch 'puzzles' into workingworking
Change-Id: I68d0d0820a30054ee67a71f36bddb6b55615b2e6
Diffstat (limited to '')
-rw-r--r--apps/plugins/puzzles/cube.c1
-rw-r--r--apps/plugins/puzzles/puzzles.h1
-rw-r--r--apps/plugins/puzzles/puzzles.make4
-rw-r--r--apps/plugins/puzzles/rbcompat.h18
-rw-r--r--apps/plugins/puzzles/rbwrappers.c2032
-rw-r--r--apps/plugins/puzzles/rockbox.c285
6 files changed, 1311 insertions, 1030 deletions
diff --git a/apps/plugins/puzzles/cube.c b/apps/plugins/puzzles/cube.c
index d0d9525..5a09648 100644
--- a/apps/plugins/puzzles/cube.c
+++ b/apps/plugins/puzzles/cube.c
@@ -361,6 +361,7 @@ static void enum_grid_squares(const game_params *params, egc_callback callback,
} else {
int row, rowlen, other, i, firstix = -1;
float theight = (float)(sqrt(3) / 2.0);
+ //float theight = 0.8660254037844386467;
for (row = 0; row < params->d1 + params->d2; row++) {
if (row < params->d2) {
diff --git a/apps/plugins/puzzles/puzzles.h b/apps/plugins/puzzles/puzzles.h
index 7201d5d..6a24ab6 100644
--- a/apps/plugins/puzzles/puzzles.h
+++ b/apps/plugins/puzzles/puzzles.h
@@ -23,6 +23,7 @@
#define lenof(array) ( sizeof(array) / sizeof(*(array)) )
+#undef STR
#define STR_INT(x) #x
#define STR(x) STR_INT(x)
diff --git a/apps/plugins/puzzles/puzzles.make b/apps/plugins/puzzles/puzzles.make
index 45326f4..537e9cf 100644
--- a/apps/plugins/puzzles/puzzles.make
+++ b/apps/plugins/puzzles/puzzles.make
@@ -17,13 +17,13 @@ OTHER_SRC += $(PUZZLES_SRC)
ifndef APP_TYPE
ROCKS += $(PUZZLES_OBJDIR)/puzzles.ovl
- PUZZLES_OUTLDS = $(PUZZLES_OBJDIR)/picutreflow.link
+ PUZZLES_OUTLDS = $(PUZZLES_OBJDIR)/puzzles.link
PUZZLES_OVLFLAGS = -T$(PUZZLES_OUTLDS) -Wl,--gc-sections -Wl,-Map,$(basename $@).map
else
ROCKS += $(PUZZLES_OBJDIR)/puzzles.rock
endif
-PUZZLESFLAGS = $(filter-out -O%,$(PLUGINFLAGS)) -O2
+PUZZLESFLAGS = $(filter-out -O%,$(PLUGINFLAGS)) -Os -Wno-unused-parameter -Wno-sign-compare
$(PUZZLES_OBJDIR)/puzzles.rock: $(PUZZLES_OBJ) $(TLSFLIB)
diff --git a/apps/plugins/puzzles/rbcompat.h b/apps/plugins/puzzles/rbcompat.h
index c733659..148aaef 100644
--- a/apps/plugins/puzzles/rbcompat.h
+++ b/apps/plugins/puzzles/rbcompat.h
@@ -7,26 +7,24 @@ int puts_wrapper(const char *s);
double sin_wrapper(double rads);
double cos_wrapper(double rads);
int vsprintf_wrapper(char *s, const char *fmt, va_list ap);
-double fabs_wrapper(double n);
-double floor_wrapper(double n);
+float fabs_wrapper(float n);
+float floor_wrapper(float n);
-double atan_wrapper(double x);
-double atan2_wrapper(double y, double x);
-double sqrt_wrapper(double x);
+float atan_wrapper(float x);
+float atan2_wrapper(float y, float x);
+float sqrt_wrapper(float x);
long strtol_wrapper(const char *nptr, char **endptr, int base);
int64_t strtoq_wrapper(const char *nptr, char **endptr, int base);
uint64_t strtouq_wrapper(const char *nptr, char **endptr, int base);
-double pow_wrapper(double x, double y);
-double scalbn_wrapper (double x, int n);
-double ceil_wrapper(double x);
+float pow_wrapper(float x, float y);
+float ceil_wrapper(float x);
size_t strspn_wrapper(const char *s1, const char *s2);
size_t strcspn_wrapper(const char *s1, const char *s2);
int sscanf_wrapper(const char *ibuf, const char *fmt, ...);
-double atof_wrapper(char *s);
+double atof_wrapper(const char *s);
double acos_wrapper(double x);
-
#define acos acos_wrapper
#define atan atan_wrapper
#define atan2 atan2_wrapper
diff --git a/apps/plugins/puzzles/rbwrappers.c b/apps/plugins/puzzles/rbwrappers.c
index c0578e1..8eca190 100644
--- a/apps/plugins/puzzles/rbwrappers.c
+++ b/apps/plugins/puzzles/rbwrappers.c
@@ -21,6 +21,7 @@ int puts_wrapper(const char *s)
return 0;
}
+/* fixed-point wrappers */
double sin_wrapper(double rads)
{
int degs = rads * 180/PI;
@@ -38,65 +39,65 @@ double cos_wrapper(double rads)
int vsprintf_wrapper(char *s, const char *fmt, va_list ap)
{
return rb->vsnprintf(s, 9999, fmt, ap);
+
}
-double fabs_wrapper(double n)
+/* Absolute value, simple calculus */
+float fabs_wrapper(float x)
{
- if(n < 0)
- return -n;
- else
- return n;
+ return (x < 0.0f) ? -x : x;
}
-double floor_wrapper(double n)
+float floor_wrapper(float n)
{
- if(n < 0)
- return ((int)n) - 1;
+ if(n < 0.0f)
+ return ((int)n - 1);
else
return (int)n;
}
-double atan_wrapper(double x);
+/* Natural logarithm.
+ Taken from glibc-2.8 */
+static const float
+ln2_hi = 6.9313812256e-01, /* 0x3f317180 */
+ln2_lo = 9.0580006145e-06, /* 0x3717f7d1 */
+two25 = 3.355443200e+07, /* 0x4c000000 */
+Lg1 = 6.6666668653e-01, /* 3F2AAAAB */
+Lg2 = 4.0000000596e-01, /* 3ECCCCCD */
+Lg3 = 2.8571429849e-01, /* 3E924925 */
+Lg4 = 2.2222198546e-01, /* 3E638E29 */
+Lg5 = 1.8183572590e-01, /* 3E3A3325 */
+Lg6 = 1.5313838422e-01, /* 3E1CD04F */
+Lg7 = 1.4798198640e-01; /* 3E178897 */
+
+static const float zero = 0.0;
+
+/* A union which permits us to convert between a float and a 32 bit
+ int. */
+
+typedef union
+{
+ float value;
+ uint32_t word;
+} ieee_float_shape_type;
-/* @(#)e_atan2.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.
- * ====================================================
- *
- */
+/* Get a 32 bit int from a float. */
-/* __ieee754_atan2(y,x)
- * Method :
- * 1. Reduce y to positive by atan2(y,x)=-atan2(-y,x).
- * 2. Reduce x to positive by (if x and y are unexceptional):
- * ARG (x+iy) = arctan(y/x) ... if x > 0,
- * ARG (x+iy) = pi - arctan[y/(-x)] ... if x < 0,
- *
- * Special cases:
- *
- * ATAN2((anything), NaN ) is NaN;
- * ATAN2(NAN , (anything) ) is NaN;
- * ATAN2(+-0, +(anything but NaN)) is +-0 ;
- * ATAN2(+-0, -(anything but NaN)) is +-pi ;
- * ATAN2(+-(anything but 0 and NaN), 0) is +-pi/2;
- * ATAN2(+-(anything but INF and NaN), +INF) is +-0 ;
- * ATAN2(+-(anything but INF and NaN), -INF) is +-pi;
- * ATAN2(+-INF,+INF ) is +-pi/4 ;
- * ATAN2(+-INF,-INF ) is +-3pi/4;
- * ATAN2(+-INF, (anything but,0,NaN, and INF)) is +-pi/2;
- *
- * Constants:
- * The hexadecimal values are the intended ones for the following
- * constants. The decimal values may be used, provided that the
- * compiler will convert from decimal to binary accurately enough
- * to produce the hexadecimal values shown.
- */
+#define GET_FLOAT_WORD(i,d) \
+do { \
+ ieee_float_shape_type gf_u; \
+ gf_u.value = (d); \
+ (i) = gf_u.word; \
+} while (0)
+
+/* Set a float from a 32 bit int. */
+
+#define SET_FLOAT_WORD(d,i) \
+do { \
+ ieee_float_shape_type sf_u; \
+ sf_u.word = (i); \
+ (d) = sf_u.value; \
+} while (0)
#ifdef ROCKBOX_LITTLE_ENDIAN
#define __HI(x) *(1+(int*)&x)
@@ -110,194 +111,752 @@ double atan_wrapper(double x);
#define __LOp(x) *(1+(int*)x)
#endif
-#ifdef __STDC__
-static const double
-#else
-static double
-#endif
-tiny = 1.0e-300,
-zero = 0.0,
-pi_o_4 = 7.8539816339744827900E-01, /* 0x3FE921FB, 0x54442D18 */
-pi_o_2 = 1.5707963267948965580E+00, /* 0x3FF921FB, 0x54442D18 */
-pi = 3.1415926535897931160E+00, /* 0x400921FB, 0x54442D18 */
-pi_lo = 1.2246467991473531772E-16; /* 0x3CA1A626, 0x33145C07 */
-
-double atan2_wrapper(double y, double x)
+static float rb_log(float x)
{
- double z;
- int k,m,hx,hy,ix,iy;
- unsigned lx,ly;
-
- hx = __HI(x); ix = hx&0x7fffffff;
- lx = __LO(x);
- hy = __HI(y); iy = hy&0x7fffffff;
- ly = __LO(y);
- if(((ix|((lx|-lx)>>31))>0x7ff00000)||
- ((iy|((ly|-ly)>>31))>0x7ff00000)) /* x or y is NaN */
- return x+y;
- if((hx-0x3ff00000|lx)==0) return atan_wrapper(y); /* x=1.0 */
- m = ((hy>>31)&1)|((hx>>30)&2); /* 2*sign(x)+sign(y) */
-
- /* when y = 0 */
- if((iy|ly)==0) {
- switch(m) {
- case 0:
- case 1: return y; /* atan(+-0,+anything)=+-0 */
- case 2: return pi+tiny;/* atan(+0,-anything) = pi */
- case 3: return -pi-tiny;/* atan(-0,-anything) =-pi */
+ float hfsq, f, s, z, R, w, t1, t2, dk;
+ int32_t k, ix, i, j;
+
+ GET_FLOAT_WORD(ix,x);
+
+ k=0;
+ if (ix < 0x00800000) { /* x < 2**-126 */
+ if ((ix&0x7fffffff)==0)
+ return -two25/(x-x); /* log(+-0)=-inf */
+ if (ix<0) return (x-x)/(x-x); /* log(-#) = NaN */
+ k -= 25; x *= two25; /* subnormal number, scale up x */
+ GET_FLOAT_WORD(ix,x);
+ }
+ if (ix >= 0x7f800000) return x+x;
+ k += (ix>>23)-127;
+ ix &= 0x007fffff;
+ i = (ix+(0x95f64<<3))&0x800000;
+ SET_FLOAT_WORD(x,ix|(i^0x3f800000)); /* normalize x or x/2 */
+ k += (i>>23);
+ f = x-(float)1.0;
+ if((0x007fffff&(15+ix))<16) { /* |f| < 2**-20 */
+ if(f==zero) {
+ if(k==0)
+ return zero;
+ else
+ {
+ dk=(float)k;
+ return dk*ln2_hi+dk*ln2_lo;
}
}
- /* when x = 0 */
- if((ix|lx)==0) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
-
- /* when x is INF */
- if(ix==0x7ff00000) {
- if(iy==0x7ff00000) {
- switch(m) {
- case 0: return pi_o_4+tiny;/* atan(+INF,+INF) */
- case 1: return -pi_o_4-tiny;/* atan(-INF,+INF) */
- case 2: return 3.0*pi_o_4+tiny;/*atan(+INF,-INF)*/
- case 3: return -3.0*pi_o_4-tiny;/*atan(-INF,-INF)*/
- }
- } else {
- switch(m) {
- case 0: return zero ; /* atan(+...,+INF) */
- case 1: return -zero ; /* atan(-...,+INF) */
- case 2: return pi+tiny ; /* atan(+...,-INF) */
- case 3: return -pi-tiny ; /* atan(-...,-INF) */
- }
- }
+ R = f*f*((float)0.5-(float)0.33333333333333333*f);
+ if(k==0)
+ return f-R;
+ else
+ {
+ dk=(float)k;
+ return dk*ln2_hi-((R-dk*ln2_lo)-f);
}
- /* when y is INF */
- if(iy==0x7ff00000) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
+ }
+ s = f/((float)2.0+f);
+ dk = (float)k;
+ z = s*s;
+ i = ix-(0x6147a<<3);
+ w = z*z;
+ j = (0x6b851<<3)-ix;
+ t1= w*(Lg2+w*(Lg4+w*Lg6));
+ t2= z*(Lg1+w*(Lg3+w*(Lg5+w*Lg7)));
+ i |= j;
+ R = t2+t1;
+ if(i>0) {
+ hfsq=(float)0.5*f*f;
+ if(k==0)
+ return f-(hfsq-s*(hfsq+R));
+ else
+ return dk*ln2_hi-((hfsq-(s*(hfsq+R)+dk*ln2_lo))-f);
+ } else {
+ if(k==0)
+ return f-s*(f-R);
+ else
+ return dk*ln2_hi-((s*(f-R)-dk*ln2_lo)-f);
+ }
+}
- /* compute y/x */
- k = (iy-ix)>>20;
- if(k > 60) z=pi_o_2+0.5*pi_lo; /* |y/x| > 2**60 */
- else if(hx<0&&k<-60) z=0.0; /* |y|/x < -2**60 */
- else z=atan_wrapper(fabs_wrapper(y/x)); /* safe to do y/x */
- switch (m) {
- case 0: return z ; /* atan(+,+) */
- case 1: __HI(z) ^= 0x80000000;
- return z ; /* atan(-,+) */
- case 2: return pi-(z-pi_lo);/* atan(+,-) */
- default: /* case 3 */
- return (z-pi_lo)-pi;/* atan(-,-) */
- }
+union ieee754_double
+ {
+ double d;
+
+ /* This is the IEEE 754 double-precision format. */
+ struct
+ {
+#ifdef ROCKBOX_BIG_ENDIAN
+ unsigned int negative:1;
+ unsigned int exponent:11;
+ /* Together these comprise the mantissa. */
+ unsigned int mantissa0:20;
+ unsigned int mantissa1:32;
+#else /* ROCKBOX_LITTLE_ENDIAN */
+ /* Together these comprise the mantissa. */
+ unsigned int mantissa1:32;
+ unsigned int mantissa0:20;
+ unsigned int exponent:11;
+ unsigned int negative:1;
+#endif /* ROCKBOX_LITTLE_ENDIAN */
+ } ieee;
+
+ /* This format makes it easier to see if a NaN is a signalling NaN. */
+ struct
+ {
+#ifdef ROCKBOX_BIG_ENDIAN
+ unsigned int negative:1;
+ unsigned int exponent:11;
+ unsigned int quiet_nan:1;
+ /* Together these comprise the mantissa. */
+ unsigned int mantissa0:19;
+ unsigned int mantissa1:32;
+#else /* ROCKBOX_LITTLE_ENDIAN */
+ /* Together these comprise the mantissa. */
+ unsigned int mantissa1:32;
+ unsigned int mantissa0:19;
+ unsigned int quiet_nan:1;
+ unsigned int exponent:11;
+ unsigned int negative:1;
+#endif /* ROCKBOX_LITTLE_ENDIAN */
+ } ieee_nan;
+ };
+
+static const volatile float TWOM100 = 7.88860905e-31;
+static const volatile float TWO127 = 1.7014118346e+38;
+
+/* Exponential function,
+ taken from glibc-2.8
+ As it uses double values and udefines some symbols,
+ it was moved to the end of the source code */
+
+#define W52 (2.22044605e-16)
+#define W55 (2.77555756e-17)
+#define W58 (3.46944695e-18)
+#define W59 (1.73472348e-18)
+#define W60 (8.67361738e-19)
+const float __exp_deltatable[178] = {
+ 0*W60, 16558714*W60, -10672149*W59, 1441652*W60,
+ -15787963*W55, 462888*W60, 7291806*W60, 1698880*W60,
+ -14375103*W58, -2021016*W60, 728829*W60, -3759654*W60,
+ 3202123*W60, -10916019*W58, -251570*W60, -1043086*W60,
+ 8207536*W60, -409964*W60, -5993931*W60, -475500*W60,
+ 2237522*W60, 324170*W60, -244117*W60, 32077*W60,
+ 123907*W60, -1019734*W60, -143*W60, 813077*W60,
+ 743345*W60, 462461*W60, 629794*W60, 2125066*W60,
+ -2339121*W60, -337951*W60, 9922067*W60, -648704*W60,
+ 149407*W60, -2687209*W60, -631608*W60, 2128280*W60,
+ -4882082*W60, 2001360*W60, 175074*W60, 2923216*W60,
+ -538947*W60, -1212193*W60, -1920926*W60, -1080577*W60,
+ 3690196*W60, 2643367*W60, 2911937*W60, 671455*W60,
+ -1128674*W60, 593282*W60, -5219347*W60, -1941490*W60,
+ 11007953*W60, 239609*W60, -2969658*W60, -1183650*W60,
+ 942998*W60, 699063*W60, 450569*W60, -329250*W60,
+ -7257875*W60, -312436*W60, 51626*W60, 555877*W60,
+ -641761*W60, 1565666*W60, 884327*W60, -10960035*W60,
+ -2004679*W60, -995793*W60, -2229051*W60, -146179*W60,
+ -510327*W60, 1453482*W60, -3778852*W60, -2238056*W60,
+ -4895983*W60, 3398883*W60, -252738*W60, 1230155*W60,
+ 346918*W60, 1109352*W60, 268941*W60, -2930483*W60,
+ -1036263*W60, -1159280*W60, 1328176*W60, 2937642*W60,
+ -9371420*W60, -6902650*W60, -1419134*W60, 1442904*W60,
+ -1319056*W60, -16369*W60, 696555*W60, -279987*W60,
+ -7919763*W60, 252741*W60, 459711*W60, -1709645*W60,
+ 354913*W60, 6025867*W60, -421460*W60, -853103*W60,
+ -338649*W60, 962151*W60, 955965*W60, 784419*W60,
+ -3633653*W60, 2277133*W60, -8847927*W52, 1223028*W60,
+ 5907079*W60, 623167*W60, 5142888*W60, 2599099*W60,
+ 1214280*W60, 4870359*W60, 593349*W60, -57705*W60,
+ 7761209*W60, -5564097*W60, 2051261*W60, 6216869*W60,
+ 4692163*W60, 601691*W60, -5264906*W60, 1077872*W60,
+ -3205949*W60, 1833082*W60, 2081746*W60, -987363*W60,
+ -1049535*W60, 2015244*W60, 874230*W60, 2168259*W60,
+ -1740124*W60, -10068269*W60, -18242*W60, -3013583*W60,
+ 580601*W60, -2547161*W60, -535689*W60, 2220815*W60,
+ 1285067*W60, 2806933*W60, -983086*W60, -1729097*W60,
+ -1162985*W60, -2561904*W60, 801988*W60, 244351*W60,
+ 1441893*W60, -7517981*W60, 271781*W60, -15021588*W60,
+ -2341588*W60, -919198*W60, 1642232*W60, 4771771*W60,
+ -1220099*W60, -3062372*W60, 628624*W60, 1278114*W60,
+ 13083513*W60, -10521925*W60, 3180310*W60, -1659307*W60,
+ 3543773*W60, 2501203*W60, 4151*W60, -340748*W60,
+ -2285625*W60, 2495202*W60
+};
+
+const double __exp_atable[355] /* __attribute__((mode(DF))) */ = {
+ 0.707722561055888932371, /* 0x0.b52d4e46605c27ffd */
+ 0.709106182438804188967, /* 0x0.b587fb96f75097ffb */
+ 0.710492508843861281234, /* 0x0.b5e2d649899167ffd */
+ 0.711881545564593931623, /* 0x0.b63dde74d36bdfffe */
+ 0.713273297897442870573, /* 0x0.b699142f945f87ffc */
+ 0.714667771153751463236, /* 0x0.b6f477909c4ea0001 */
+ 0.716064970655995725059, /* 0x0.b75008aec758f8004 */
+ 0.717464901723956938193, /* 0x0.b7abc7a0eea7e0002 */
+ 0.718867569715736398602, /* 0x0.b807b47e1586c7ff8 */
+ 0.720272979947266023271, /* 0x0.b863cf5d10e380003 */
+ 0.721681137825144314297, /* 0x0.b8c01855195c37ffb */
+ 0.723092048691992950199, /* 0x0.b91c8f7d213740004 */
+ 0.724505717938892290800, /* 0x0.b97934ec5002d0007 */
+ 0.725922150953176470431, /* 0x0.b9d608b9c92ea7ffc */
+ 0.727341353138962865022, /* 0x0.ba330afcc29e98003 */
+ 0.728763329918453162104, /* 0x0.ba903bcc8618b7ffc */
+ 0.730188086709957051568, /* 0x0.baed9b40591ba0000 */
+ 0.731615628948127705309, /* 0x0.bb4b296f931e30002 */
+ 0.733045962086486091436, /* 0x0.bba8e671a05617ff9 */
+ 0.734479091556371366251, /* 0x0.bc06d25dd49568001 */
+ 0.735915022857225542529, /* 0x0.bc64ed4bce8f6fff9 */
+ 0.737353761441304711410, /* 0x0.bcc33752f915d7ff9 */
+ 0.738795312814142124419, /* 0x0.bd21b08af98e78005 */
+ 0.740239682467211168593, /* 0x0.bd80590b65e9a8000 */
+ 0.741686875913991849885, /* 0x0.bddf30ebec4a10000 */
+ 0.743136898669507939299, /* 0x0.be3e38443c84e0007 */
+ 0.744589756269486091620, /* 0x0.be9d6f2c1d32a0002 */
+ 0.746045454254026796384, /* 0x0.befcd5bb59baf8004 */
+ 0.747503998175051087583, /* 0x0.bf5c6c09ca84c0003 */
+ 0.748965393601880857739, /* 0x0.bfbc322f5b18b7ff8 */
+ 0.750429646104262104698, /* 0x0.c01c2843f776fffff */
+ 0.751896761271877989160, /* 0x0.c07c4e5fa18b88002 */
+ 0.753366744698445112140, /* 0x0.c0dca49a5fb18fffd */
+ 0.754839601988627206827, /* 0x0.c13d2b0c444db0005 */
+ 0.756315338768691947122, /* 0x0.c19de1cd798578006 */
+ 0.757793960659406629066, /* 0x0.c1fec8f623723fffd */
+ 0.759275473314173443536, /* 0x0.c25fe09e8a0f47ff8 */
+ 0.760759882363831851927, /* 0x0.c2c128dedc88f8000 */
+ 0.762247193485956486805, /* 0x0.c322a1cf7d6e7fffa */
+ 0.763737412354726363781, /* 0x0.c3844b88cb9347ffc */
+ 0.765230544649828092739, /* 0x0.c3e626232bd8f7ffc */
+ 0.766726596071518051729, /* 0x0.c44831b719bf18002 */
+ 0.768225572321911687194, /* 0x0.c4aa6e5d12d078001 */
+ 0.769727479119219348810, /* 0x0.c50cdc2da64a37ffb */
+ 0.771232322196981678892, /* 0x0.c56f7b41744490001 */
+ 0.772740107296721268087, /* 0x0.c5d24bb1259e70004 */
+ 0.774250840160724651565, /* 0x0.c6354d95640dd0007 */
+ 0.775764526565368872643, /* 0x0.c6988106fec447fff */
+ 0.777281172269557396602, /* 0x0.c6fbe61eb1bd0ffff */
+ 0.778800783068235302750, /* 0x0.c75f7cf560942fffc */
+ 0.780323364758801041312, /* 0x0.c7c345a3f1983fffe */
+ 0.781848923151573727006, /* 0x0.c8274043594cb0002 */
+ 0.783377464064598849602, /* 0x0.c88b6cec94b3b7ff9 */
+ 0.784908993312207869935, /* 0x0.c8efcbb89cba27ffe */
+ 0.786443516765346961618, /* 0x0.c9545cc0a88c70003 */
+ 0.787981040257604625744, /* 0x0.c9b9201dc643bfffa */
+ 0.789521569657452682047, /* 0x0.ca1e15e92a5410007 */
+ 0.791065110849462849192, /* 0x0.ca833e3c1ae510005 */
+ 0.792611669712891875319, /* 0x0.cae8992fd84667ffd */
+ 0.794161252150049179450, /* 0x0.cb4e26ddbc207fff8 */
+ 0.795713864077794763584, /* 0x0.cbb3e75f301b60003 */
+ 0.797269511407239561694, /* 0x0.cc19dacd978cd8002 */
+ 0.798828200086368567220, /* 0x0.cc8001427e55d7ffb */
+ 0.800389937624300440456, /* 0x0.cce65ade24d360006 */
+ 0.801954725261124767840, /* 0x0.cd4ce7a5de839fffb */
+ 0.803522573691593189330, /* 0x0.cdb3a7c79a678fffd */
+ 0.805093487311204114563, /* 0x0.ce1a9b563965ffffc */
+ 0.806667472122675088819, /* 0x0.ce81c26b838db8000 */
+ 0.808244534127439906441, /* 0x0.cee91d213f8428002 */
+ 0.809824679342317166307, /* 0x0.cf50ab9144d92fff9 */
+ 0.811407913793616542005, /* 0x0.cfb86dd5758c2ffff */
+ 0.812994243520784198882, /* 0x0.d0206407c20e20005 */
+ 0.814583674571603966162, /* 0x0.d0888e4223facfff9 */
+ 0.816176213022088536960, /* 0x0.d0f0ec9eb3f7c8002 */
+ 0.817771864936188586101, /* 0x0.d1597f377d6768002 */
+ 0.819370636400374108252, /* 0x0.d1c24626a46eafff8 */
+ 0.820972533518165570298, /* 0x0.d22b41865ff1e7ff9 */
+ 0.822577562404315121269, /* 0x0.d2947170f32ec7ff9 */
+ 0.824185729164559344159, /* 0x0.d2fdd60097795fff8 */
+ 0.825797039949601741075, /* 0x0.d3676f4fb796d0001 */
+ 0.827411500902565544264, /* 0x0.d3d13d78b5f68fffb */
+ 0.829029118181348834154, /* 0x0.d43b40960546d8001 */
+ 0.830649897953322891022, /* 0x0.d4a578c222a058000 */
+ 0.832273846408250750368, /* 0x0.d50fe617a3ba78005 */
+ 0.833900969738858188772, /* 0x0.d57a88b1218e90002 */
+ 0.835531274148056613016, /* 0x0.d5e560a94048f8006 */
+ 0.837164765846411529371, /* 0x0.d6506e1aac8078003 */
+ 0.838801451086016225394, /* 0x0.d6bbb1204074e0001 */
+ 0.840441336100884561780, /* 0x0.d72729d4c28518004 */
+ 0.842084427144139224814, /* 0x0.d792d8530e12b0001 */
+ 0.843730730487052604790, /* 0x0.d7febcb61273e7fff */
+ 0.845380252404570153833, /* 0x0.d86ad718c308dfff9 */
+ 0.847032999194574087728, /* 0x0.d8d727962c69d7fff */
+ 0.848688977161248581090, /* 0x0.d943ae49621ce7ffb */
+ 0.850348192619261200615, /* 0x0.d9b06b4d832ef8005 */
+ 0.852010651900976245816, /* 0x0.da1d5ebdc22220005 */
+ 0.853676361342631029337, /* 0x0.da8a88b555baa0006 */
+ 0.855345327311054837175, /* 0x0.daf7e94f965f98004 */
+ 0.857017556155879489641, /* 0x0.db6580a7c98f7fff8 */
+ 0.858693054267390953857, /* 0x0.dbd34ed9617befff8 */
+ 0.860371828028939855647, /* 0x0.dc4153ffc8b65fff9 */
+ 0.862053883854957292436, /* 0x0.dcaf90368bfca8004 */
+ 0.863739228154875360306, /* 0x0.dd1e0399328d87ffe */
+ 0.865427867361348468455, /* 0x0.dd8cae435d303fff9 */
+ 0.867119807911702289458, /* 0x0.ddfb9050b1cee8006 */
+ 0.868815056264353846599, /* 0x0.de6aa9dced8448001 */
+ 0.870513618890481399881, /* 0x0.ded9fb03db7320006 */
+ 0.872215502247877139094, /* 0x0.df4983e1380657ff8 */
+ 0.873920712852848668986, /* 0x0.dfb94490ffff77ffd */
+ 0.875629257204025623884, /* 0x0.e0293d2f1cb01fff9 */
+ 0.877341141814212965880, /* 0x0.e0996dd786fff0007 */
+ 0.879056373217612985183, /* 0x0.e109d6a64f5d57ffc */
+ 0.880774957955916648615, /* 0x0.e17a77b78e72a7ffe */
+ 0.882496902590150900078, /* 0x0.e1eb5127722cc7ff8 */
+ 0.884222213673356738383, /* 0x0.e25c63121fb0c8006 */
+ 0.885950897802399772740, /* 0x0.e2cdad93ec5340003 */
+ 0.887682961567391237685, /* 0x0.e33f30c925fb97ffb */
+ 0.889418411575228162725, /* 0x0.e3b0ecce2d05ffff9 */
+ 0.891157254447957902797, /* 0x0.e422e1bf727718006 */
+ 0.892899496816652704641, /* 0x0.e4950fb9713fc7ffe */
+ 0.894645145323828439008, /* 0x0.e50776d8b0e60fff8 */
+ 0.896394206626591749641, /* 0x0.e57a1739c8fadfffc */
+ 0.898146687421414902124, /* 0x0.e5ecf0f97c5798007 */
+ 0.899902594367530173098, /* 0x0.e660043464e378005 */
+ 0.901661934163603406867, /* 0x0.e6d3510747e150006 */
+ 0.903424713533971135418, /* 0x0.e746d78f06cd97ffd */
+ 0.905190939194458810123, /* 0x0.e7ba97e879c91fffc */
+ 0.906960617885092856864, /* 0x0.e82e92309390b0007 */
+ 0.908733756358986566306, /* 0x0.e8a2c6845544afffa */
+ 0.910510361377119825629, /* 0x0.e9173500c8abc7ff8 */
+ 0.912290439722343249336, /* 0x0.e98bddc30f98b0002 */
+ 0.914073998177417412765, /* 0x0.ea00c0e84bc4c7fff */
+ 0.915861043547953501680, /* 0x0.ea75de8db8094fffe */
+ 0.917651582652244779397, /* 0x0.eaeb36d09d3137ffe */
+ 0.919445622318405764159, /* 0x0.eb60c9ce4ed3dffff */
+ 0.921243169397334638073, /* 0x0.ebd697a43995b0007 */
+ 0.923044230737526172328, /* 0x0.ec4ca06fc7768fffa */
+ 0.924848813220121135342, /* 0x0.ecc2e44e865b6fffb */
+ 0.926656923710931002014, /* 0x0.ed39635df34e70006 */
+ 0.928468569126343790092, /* 0x0.edb01dbbc2f5b7ffa */
+ 0.930283756368834757725, /* 0x0.ee2713859aab57ffa */
+ 0.932102492359406786818, /* 0x0.ee9e44d9342870004 */
+ 0.933924784042873379360, /* 0x0.ef15b1d4635438005 */
+ 0.935750638358567643520, /* 0x0.ef8d5a94f60f50007 */
+ 0.937580062297704630580, /* 0x0.f0053f38f345cffff */
+ 0.939413062815381727516, /* 0x0.f07d5fde3a2d98001 */
+ 0.941249646905368053689, /* 0x0.f0f5bca2d481a8004 */
+ 0.943089821583810716806, /* 0x0.f16e55a4e497d7ffe */
+ 0.944933593864477061592, /* 0x0.f1e72b028a2827ffb */
+ 0.946780970781518460559, /* 0x0.f2603cd9fb5430001 */
+ 0.948631959382661205081, /* 0x0.f2d98b497d2a87ff9 */
+ 0.950486566729423554277, /* 0x0.f353166f63e3dffff */
+ 0.952344799896018723290, /* 0x0.f3ccde6a11ae37ffe */
+ 0.954206665969085765512, /* 0x0.f446e357f66120000 */
+ 0.956072172053890279009, /* 0x0.f4c12557964f0fff9 */
+ 0.957941325265908139014, /* 0x0.f53ba48781046fffb */
+ 0.959814132734539637840, /* 0x0.f5b66106555d07ffa */
+ 0.961690601603558903308, /* 0x0.f6315af2c2027fffc */
+ 0.963570739036113010927, /* 0x0.f6ac926b8aeb80004 */
+ 0.965454552202857141381, /* 0x0.f728078f7c5008002 */
+ 0.967342048278315158608, /* 0x0.f7a3ba7d66a908001 */
+ 0.969233234469444204768, /* 0x0.f81fab543e1897ffb */
+ 0.971128118008140250896, /* 0x0.f89bda33122c78007 */
+ 0.973026706099345495256, /* 0x0.f9184738d4cf97ff8 */
+ 0.974929006031422851235, /* 0x0.f994f284d3a5c0008 */
+ 0.976835024947348973265, /* 0x0.fa11dc35bc7820002 */
+ 0.978744770239899142285, /* 0x0.fa8f046b4fb7f8007 */
+ 0.980658249138918636210, /* 0x0.fb0c6b449ab1cfff9 */
+ 0.982575468959622777535, /* 0x0.fb8a10e1088fb7ffa */
+ 0.984496437054508843888, /* 0x0.fc07f5602d79afffc */
+ 0.986421160608523028820, /* 0x0.fc8618e0e55e47ffb */
+ 0.988349647107594098099, /* 0x0.fd047b83571b1fffa */
+ 0.990281903873210800357, /* 0x0.fd831d66f4c018002 */
+ 0.992217938695037382475, /* 0x0.fe01fead3320bfff8 */
+ 0.994157757657894713987, /* 0x0.fe811f703491e8006 */
+ 0.996101369488558541238, /* 0x0.ff007fd5744490005 */
+ 0.998048781093141101932, /* 0x0.ff801ffa9b9280007 */
+ 1.000000000000000000000, /* 0x1.00000000000000000 */
+ 1.001955033605393285965, /* 0x1.0080200565d29ffff */
+ 1.003913889319761887310, /* 0x1.0100802aa0e80fff0 */
+ 1.005876574715736104818, /* 0x1.01812090377240007 */
+ 1.007843096764807100351, /* 0x1.020201541aad7fff6 */
+ 1.009813464316352327214, /* 0x1.0283229c4c9820007 */
+ 1.011787683565730677817, /* 0x1.030484836910a000e */
+ 1.013765762469146736174, /* 0x1.0386272b9c077fffe */
+ 1.015747708536026694351, /* 0x1.04080ab526304fff0 */
+ 1.017733529475172815584, /* 0x1.048a2f412375ffff0 */
+ 1.019723232714418781378, /* 0x1.050c94ef7ad5e000a */
+ 1.021716825883923762690, /* 0x1.058f3be0f1c2d0004 */
+ 1.023714316605201180057, /* 0x1.06122436442e2000e */
+ 1.025715712440059545995, /* 0x1.06954e0fec63afff2 */
+ 1.027721021151397406936, /* 0x1.0718b98f41c92fff6 */
+ 1.029730250269221158939, /* 0x1.079c66d49bb2ffff1 */
+ 1.031743407506447551857, /* 0x1.082056011a9230009 */
+ 1.033760500517691527387, /* 0x1.08a487359ebd50002 */
+ 1.035781537016238873464, /* 0x1.0928fa93490d4fff3 */
+ 1.037806524719013578963, /* 0x1.09adb03b3e5b3000d */
+ 1.039835471338248051878, /* 0x1.0a32a84e9e5760004 */
+ 1.041868384612101516848, /* 0x1.0ab7e2eea5340ffff */
+ 1.043905272300907460835, /* 0x1.0b3d603ca784f0009 */
+ 1.045946142174331239262, /* 0x1.0bc3205a042060000 */
+ 1.047991002016745332165, /* 0x1.0c4923682a086fffe */
+ 1.050039859627715177527, /* 0x1.0ccf698898f3a000d */
+ 1.052092722826109660856, /* 0x1.0d55f2dce5d1dfffb */
+ 1.054149599440827866881, /* 0x1.0ddcbf86b09a5fff6 */
+ 1.056210497317612961855, /* 0x1.0e63cfa7abc97fffd */
+ 1.058275424318780855142, /* 0x1.0eeb23619c146fffb */
+ 1.060344388322010722446, /* 0x1.0f72bad65714bffff */
+ 1.062417397220589476718, /* 0x1.0ffa9627c38d30004 */
+ 1.064494458915699715017, /* 0x1.1082b577d0eef0003 */
+ 1.066575581342167566880, /* 0x1.110b18e893a90000a */
+ 1.068660772440545025953, /* 0x1.1193c09c267610006 */
+ 1.070750040138235936705, /* 0x1.121cacb4959befff6 */
+ 1.072843392435016474095, /* 0x1.12a5dd543cf36ffff */
+ 1.074940837302467588937, /* 0x1.132f529d59552000b */
+ 1.077042382749654914030, /* 0x1.13b90cb250d08fff5 */
+ 1.079148036789447484528, /* 0x1.14430bb58da3dfff9 */
+ 1.081257807444460983297, /* 0x1.14cd4fc984c4a000e */
+ 1.083371702785017154417, /* 0x1.1557d910df9c7000e */
+ 1.085489730853784307038, /* 0x1.15e2a7ae292d30002 */
+ 1.087611899742884524772, /* 0x1.166dbbc422d8c0004 */
+ 1.089738217537583819804, /* 0x1.16f9157586772ffff */
+ 1.091868692357631731528, /* 0x1.1784b4e533cacfff0 */
+ 1.094003332327482702577, /* 0x1.18109a360fc23fff2 */
+ 1.096142145591650907149, /* 0x1.189cc58b155a70008 */
+ 1.098285140311341168136, /* 0x1.1929370751ea50002 */
+ 1.100432324652149906842, /* 0x1.19b5eecdd79cefff0 */
+ 1.102583706811727015711, /* 0x1.1a42ed01dbdba000e */
+ 1.104739294993289488947, /* 0x1.1ad031c69a2eafff0 */
+ 1.106899097422573863281, /* 0x1.1b5dbd3f66e120003 */
+ 1.109063122341542140286, /* 0x1.1beb8f8fa8150000b */
+ 1.111231377994659874592, /* 0x1.1c79a8dac6ad0fff4 */
+ 1.113403872669181282605, /* 0x1.1d0809445a97ffffc */
+ 1.115580614653132185460, /* 0x1.1d96b0effc9db000e */
+ 1.117761612217810673898, /* 0x1.1e25a001332190000 */
+ 1.119946873713312474002, /* 0x1.1eb4d69bdb2a9fff1 */
+ 1.122136407473298902480, /* 0x1.1f4454e3bfae00006 */
+ 1.124330221845670330058, /* 0x1.1fd41afcbb48bfff8 */
+ 1.126528325196519908506, /* 0x1.2064290abc98c0001 */
+ 1.128730725913251964394, /* 0x1.20f47f31c9aa7000f */
+ 1.130937432396844410880, /* 0x1.21851d95f776dfff0 */
+ 1.133148453059692917203, /* 0x1.2216045b6784efffa */
+ 1.135363796355857157764, /* 0x1.22a733a6692ae0004 */
+ 1.137583470716100553249, /* 0x1.2338ab9b3221a0004 */
+ 1.139807484614418608939, /* 0x1.23ca6c5e27aadfff7 */
+ 1.142035846532929888057, /* 0x1.245c7613b7f6c0004 */
+ 1.144268564977221958089, /* 0x1.24eec8e06b035000c */
+ 1.146505648458203463465, /* 0x1.258164e8cea85fff8 */
+ 1.148747105501412235671, /* 0x1.26144a5180d380009 */
+ 1.150992944689175123667, /* 0x1.26a7793f5de2efffa */
+ 1.153243174560058870217, /* 0x1.273af1d712179000d */
+ 1.155497803703682491111, /* 0x1.27ceb43d81d42fff1 */
+ 1.157756840726344771440, /* 0x1.2862c097a3d29000c */
+ 1.160020294239811677834, /* 0x1.28f7170a74cf4fff1 */
+ 1.162288172883275239058, /* 0x1.298bb7bb0faed0004 */
+ 1.164560485298402170388, /* 0x1.2a20a2ce920dffff4 */
+ 1.166837240167474476460, /* 0x1.2ab5d86a4631ffff6 */
+ 1.169118446164539637555, /* 0x1.2b4b58b36d5220009 */
+ 1.171404112007080167155, /* 0x1.2be123cf786790002 */
+ 1.173694246390975415341, /* 0x1.2c7739e3c0aac000d */
+ 1.175988858069749065617, /* 0x1.2d0d9b15deb58fff6 */
+ 1.178287955789017793514, /* 0x1.2da4478b627040002 */
+ 1.180591548323240091978, /* 0x1.2e3b3f69fb794fffc */
+ 1.182899644456603782686, /* 0x1.2ed282d76421d0004 */
+ 1.185212252993012693694, /* 0x1.2f6a11f96c685fff3 */
+ 1.187529382762033236513, /* 0x1.3001ecf60082ffffa */
+ 1.189851042595508889847, /* 0x1.309a13f30f28a0004 */
+ 1.192177241354644978669, /* 0x1.31328716a758cfff7 */
+ 1.194507987909589896687, /* 0x1.31cb4686e1e85fffb */
+ 1.196843291137896336843, /* 0x1.32645269dfd04000a */
+ 1.199183159977805113226, /* 0x1.32fdaae604c39000f */
+ 1.201527603343041317132, /* 0x1.339750219980dfff3 */
+ 1.203876630171082595692, /* 0x1.3431424300e480007 */
+ 1.206230249419600664189, /* 0x1.34cb8170b3fee000e */
+ 1.208588470077065268869, /* 0x1.35660dd14dbd4fffc */
+ 1.210951301134513435915, /* 0x1.3600e78b6bdfc0005 */
+ 1.213318751604272271958, /* 0x1.369c0ec5c38ebfff2 */
+ 1.215690830512196507537, /* 0x1.373783a718d29000f */
+ 1.218067546930756250870, /* 0x1.37d3465662f480007 */
+ 1.220448909901335365929, /* 0x1.386f56fa770fe0008 */
+ 1.222834928513994334780, /* 0x1.390bb5ba5fc540004 */
+ 1.225225611877684750397, /* 0x1.39a862bd3c7a8fff3 */
+ 1.227620969111500981433, /* 0x1.3a455e2a37bcafffd */
+ 1.230021009336254911271, /* 0x1.3ae2a8287dfbefff6 */
+ 1.232425741726685064472, /* 0x1.3b8040df76f39fffa */
+ 1.234835175450728295084, /* 0x1.3c1e287682e48fff1 */
+ 1.237249319699482263931, /* 0x1.3cbc5f151b86bfff8 */
+ 1.239668183679933477545, /* 0x1.3d5ae4e2cc0a8000f */
+ 1.242091776620540377629, /* 0x1.3df9ba07373bf0006 */
+ 1.244520107762172811399, /* 0x1.3e98deaa0d8cafffe */
+ 1.246953186383919165383, /* 0x1.3f3852f32973efff0 */
+ 1.249391019292643401078, /* 0x1.3fd816ffc72b90001 */
+ 1.251833623164381181797, /* 0x1.40782b17863250005 */
+ 1.254280999953110153911, /* 0x1.41188f42caf400000 */
+ 1.256733161434815393410, /* 0x1.41b943b42945bfffd */
+ 1.259190116985283935980, /* 0x1.425a4893e5f10000a */
+ 1.261651875958665236542, /* 0x1.42fb9e0a2df4c0009 */
+ 1.264118447754797758244, /* 0x1.439d443f608c4fff9 */
+ 1.266589841787181258708, /* 0x1.443f3b5bebf850008 */
+ 1.269066067469190262045, /* 0x1.44e183883e561fff7 */
+ 1.271547134259576328224, /* 0x1.45841cecf7a7a0001 */
+ 1.274033051628237434048, /* 0x1.462707b2c43020009 */
+ 1.276523829025464573684, /* 0x1.46ca44023aa410007 */
+ 1.279019475999373156531, /* 0x1.476dd2045d46ffff0 */
+ 1.281520002043128991825, /* 0x1.4811b1e1f1f19000b */
+ 1.284025416692967214122, /* 0x1.48b5e3c3edd74fff4 */
+ 1.286535729509738823464, /* 0x1.495a67d3613c8fff7 */
+ 1.289050950070396384145, /* 0x1.49ff3e396e19d000b */
+ 1.291571087985403654081, /* 0x1.4aa4671f5b401fff1 */
+ 1.294096152842774794011, /* 0x1.4b49e2ae56d19000d */
+ 1.296626154297237043484, /* 0x1.4befb10fd84a3fff4 */
+ 1.299161101984141142272, /* 0x1.4c95d26d41d84fff8 */
+ 1.301701005575179204100, /* 0x1.4d3c46f01d9f0fff3 */
+ 1.304245874766450485904, /* 0x1.4de30ec21097d0003 */
+ 1.306795719266019562007, /* 0x1.4e8a2a0ccce3d0002 */
+ 1.309350548792467483458, /* 0x1.4f3198fa10346fff5 */
+ 1.311910373099227200545, /* 0x1.4fd95bb3be8cffffd */
+ 1.314475201942565174546, /* 0x1.50817263bf0e5fffb */
+ 1.317045045107389400535, /* 0x1.5129dd3418575000e */
+ 1.319619912422941299109, /* 0x1.51d29c4f01c54ffff */
+ 1.322199813675649204855, /* 0x1.527bafde83a310009 */
+ 1.324784758729532718739, /* 0x1.5325180cfb8b3fffd */
+ 1.327374757430096474625, /* 0x1.53ced504b2bd0fff4 */
+ 1.329969819671041886272, /* 0x1.5478e6f02775e0001 */
+ 1.332569955346704748651, /* 0x1.55234df9d8a59fff8 */
+ 1.335175174370685002822, /* 0x1.55ce0a4c5a6a9fff6 */
+ 1.337785486688218616860, /* 0x1.56791c1263abefff7 */
+ 1.340400902247843806217, /* 0x1.57248376aef21fffa */
+ 1.343021431036279800211, /* 0x1.57d040a420c0bfff3 */
+ 1.345647083048053138662, /* 0x1.587c53c5a630f0002 */
+ 1.348277868295411074918, /* 0x1.5928bd063fd7bfff9 */
+ 1.350913796821875845231, /* 0x1.59d57c9110ad60006 */
+ 1.353554878672557082439, /* 0x1.5a8292913d68cfffc */
+ 1.356201123929036356254, /* 0x1.5b2fff3212db00007 */
+ 1.358852542671913132777, /* 0x1.5bddc29edcc06fff3 */
+ 1.361509145047255398051, /* 0x1.5c8bdd032ed16000f */
+ 1.364170941142184734180, /* 0x1.5d3a4e8a5bf61fff4 */
+ 1.366837941171020309735, /* 0x1.5de9176042f1effff */
+ 1.369510155261156381121, /* 0x1.5e9837b062f4e0005 */
+ 1.372187593620959988833, /* 0x1.5f47afa69436cfff1 */
+ 1.374870266463378287715, /* 0x1.5ff77f6eb3f8cfffd */
+ 1.377558184010425845733, /* 0x1.60a7a734a9742fff9 */
+ 1.380251356531521533853, /* 0x1.6158272490016000c */
+ 1.382949794301995272203, /* 0x1.6208ff6a8978a000f */
+ 1.385653507605306700170, /* 0x1.62ba3032c0a280004 */
+ 1.388362506772382154503, /* 0x1.636bb9a994784000f */
+ 1.391076802081129493127, /* 0x1.641d9bfb29a7bfff6 */
+ 1.393796403973427855412, /* 0x1.64cfd7545928b0002 */
+ 1.396521322756352656542, /* 0x1.65826be167badfff8 */
+ 1.399251568859207761660, /* 0x1.663559cf20826000c */
+ 1.401987152677323100733, /* 0x1.66e8a14a29486fffc */
+ 1.404728084651919228815, /* 0x1.679c427f5a4b6000b */
+ 1.407474375243217723560, /* 0x1.68503d9ba0add000f */
+ 1.410226034922914983815, /* 0x1.690492cbf6303fff9 */
+ 1.412983074197955213304, /* 0x1.69b9423d7b548fff6 */
+};
+
+/* All floating-point numbers can be put in one of these categories. */
+enum
+ {
+ FP_NAN,
+# define FP_NAN FP_NAN
+ FP_INFINITE,
+# define FP_INFINITE FP_INFINITE
+ FP_ZERO,
+# define FP_ZERO FP_ZERO
+ FP_SUBNORMAL,
+# define FP_SUBNORMAL FP_SUBNORMAL
+ FP_NORMAL
+# define FP_NORMAL FP_NORMAL
+ };
+
+
+int
+__fpclassifyf (float x)
+{
+ uint32_t wx;
+ int retval = FP_NORMAL;
+
+ GET_FLOAT_WORD (wx, x);
+ wx &= 0x7fffffff;
+ if (wx == 0)
+ retval = FP_ZERO;
+ else if (wx < 0x800000)
+ retval = FP_SUBNORMAL;
+ else if (wx >= 0x7f800000)
+ retval = wx > 0x7f800000 ? FP_NAN : FP_INFINITE;
+
+ return retval;
}
-/* @(#)s_atan.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.
- * ====================================================
- *
- */
-/* atan(x)
- * Method
- * 1. Reduce x to positive by atan(x) = -atan(-x).
- * 2. According to the integer k=4t+0.25 chopped, t=x, the argument
- * is further reduced to one of the following intervals and the
- * arctangent of t is evaluated by the corresponding formula:
- *
- * [0,7/16] atan(x) = t-t^3*(a1+t^2*(a2+...(a10+t^2*a11)...)
- * [7/16,11/16] atan(x) = atan(1/2) + atan( (t-0.5)/(1+t/2) )
- * [11/16.19/16] atan(x) = atan( 1 ) + atan( (t-1)/(1+t) )
- * [19/16,39/16] atan(x) = atan(3/2) + atan( (t-1.5)/(1+1.5t) )
- * [39/16,INF] atan(x) = atan(INF) + atan( -1/t )
- *
- * Constants:
- * The hexadecimal values are the intended ones for the following
- * constants. The decimal values may be used, provided that the
- * compiler will convert from decimal to binary accurately enough
- * to produce the hexadecimal values shown.
- */
+int
+__isinff (float x)
+{
+ int32_t ix,t;
+ GET_FLOAT_WORD(ix,x);
+ t = ix & 0x7fffffff;
+ t ^= 0x7f800000;
+ t |= -t;
+ return ~(t >> 31) & (ix >> 30);
+}
-#ifdef __STDC__
-static const double atanhi[] = {
-#else
-static double atanhi[] = {
+/* Return nonzero value if arguments are unordered. */
+#define fpclassify(x) \
+ (sizeof (x) == sizeof (float) ? __fpclassifyf (x) : __fpclassifyf (x))
+
+#ifndef isunordered
+#define isunordered(u, v) \
+ (__extension__ \
+ ({ __typeof__(u) __u = (u); __typeof__(v) __v = (v); \
+ fpclassify (__u) == FP_NAN || fpclassify (__v) == FP_NAN; }))
#endif
- 4.63647609000806093515e-01, /* atan(0.5)hi 0x3FDDAC67, 0x0561BB4F */
- 7.85398163397448278999e-01, /* atan(1.0)hi 0x3FE921FB, 0x54442D18 */
- 9.82793723247329054082e-01, /* atan(1.5)hi 0x3FEF730B, 0xD281F69B */
- 1.57079632679489655800e+00, /* atan(inf)hi 0x3FF921FB, 0x54442D18 */
-};
-#ifdef __STDC__
-static const double atanlo[] = {
-#else
-static double atanlo[] = {
+/* Return nonzero value if X is less than Y. */
+#ifndef isless
+#define isless(x, y) \
+ (__extension__ \
+ ({ __typeof__(x) __x = (x); __typeof__(y) __y = (y); \
+ !isunordered (__x, __y) && __x < __y; }))
#endif
- 2.26987774529616870924e-17, /* atan(0.5)lo 0x3C7A2B7F, 0x222F65E2 */
- 3.06161699786838301793e-17, /* atan(1.0)lo 0x3C81A626, 0x33145C07 */
- 1.39033110312309984516e-17, /* atan(1.5)lo 0x3C700788, 0x7AF0CBBD */
- 6.12323399573676603587e-17, /* atan(inf)lo 0x3C91A626, 0x33145C07 */
-};
-#ifdef __STDC__
-static const double aT[] = {
-#else
-static double aT[] = {
+/* Return nonzero value if X is greater than Y. */
+#ifndef isgreater
+#define isgreater(x, y) \
+ (__extension__ \
+ ({ __typeof__(x) __x = (x); __typeof__(y) __y = (y); \
+ !isunordered (__x, __y) && __x > __y; }))
#endif
- 3.33333333333329318027e-01, /* 0x3FD55555, 0x5555550D */
- -1.99999999998764832476e-01, /* 0xBFC99999, 0x9998EBC4 */
- 1.42857142725034663711e-01, /* 0x3FC24924, 0x920083FF */
- -1.11111104054623557880e-01, /* 0xBFBC71C6, 0xFE231671 */
- 9.09088713343650656196e-02, /* 0x3FB745CD, 0xC54C206E */
- -7.69187620504482999495e-02, /* 0xBFB3B0F2, 0xAF749A6D */
- 6.66107313738753120669e-02, /* 0x3FB10D66, 0xA0D03D51 */
- -5.83357013379057348645e-02, /* 0xBFADDE2D, 0x52DEFD9A */
- 4.97687799461593236017e-02, /* 0x3FA97B4B, 0x24760DEB */
- -3.65315727442169155270e-02, /* 0xBFA2B444, 0x2C6A6C2F */
- 1.62858201153657823623e-02, /* 0x3F90AD3A, 0xE322DA11 */
-};
-#ifdef __STDC__
- static const double
-#else
- static double
+float rb_exp(float x)
+{
+ static const float himark = 88.72283935546875;
+ static const float lomark = -103.972084045410;
+ /* Check for usual case. */
+ if (isless (x, himark) && isgreater (x, lomark))
+ {
+ static const float THREEp42 = 13194139533312.0;
+ static const float THREEp22 = 12582912.0;
+ /* 1/ln(2). */
+#undef M_1_LN2
+ static const float M_1_LN2 = 1.44269502163f;
+ /* ln(2) */
+#undef M_LN2
+ static const double M_LN2 = .6931471805599452862;
+
+ int tval;
+ double x22, t, result, dx;
+ float n, delta;
+ union ieee754_double ex2_u;
+#ifndef ROCKBOX
+ fenv_t oldenv;
+
+ feholdexcept (&oldenv);
+#endif
+
+#ifdef FE_TONEAREST
+ fesetround (FE_TONEAREST);
#endif
-one = 1.0,
-huge = 1.0e300;
-double atan_wrapper(double x)
+ /* Calculate n. */
+ n = x * M_1_LN2 + THREEp22;
+ n -= THREEp22;
+ dx = x - n*M_LN2;
+
+ /* Calculate t/512. */
+ t = dx + THREEp42;
+ t -= THREEp42;
+ dx -= t;
+
+ /* Compute tval = t. */
+ tval = (int) (t * 512.0);
+
+ if (t >= 0)
+ delta = - __exp_deltatable[tval];
+ else
+ delta = __exp_deltatable[-tval];
+
+ /* Compute ex2 = 2^n e^(t/512+delta[t]). */
+ ex2_u.d = __exp_atable[tval+177];
+ ex2_u.ieee.exponent += (int) n;
+
+ /* Approximate e^(dx+delta) - 1, using a second-degree polynomial,
+ with maximum error in [-2^-10-2^-28,2^-10+2^-28]
+ less than 5e-11. */
+ x22 = (0.5000000496709180453 * dx + 1.0000001192102037084) * dx + delta;
+
+ /* Return result. */
+#ifndef ROCKBOX
+ fesetenv (&oldenv);
+#endif
+
+ result = x22 * ex2_u.d + ex2_u.d;
+ return (float) result;
+ }
+ /* Exceptional cases: */
+ else if (isless (x, himark))
+ {
+ if (__isinff (x))
+ /* e^-inf == 0, with no error. */
+ return 0;
+ else
+ /* Underflow */
+ return TWOM100 * TWOM100;
+ }
+ else
+ /* Return x, if x is a NaN or Inf; or overflow, otherwise. */
+ return TWO127*x;
+}
+
+/* Arc tangent,
+ taken from glibc-2.8. */
+
+static const float atanhi[] = {
+ 4.6364760399e-01, /* atan(0.5)hi 0x3eed6338 */
+ 7.8539812565e-01, /* atan(1.0)hi 0x3f490fda */
+ 9.8279368877e-01, /* atan(1.5)hi 0x3f7b985e */
+ 1.5707962513e+00, /* atan(inf)hi 0x3fc90fda */
+};
+
+static const float atanlo[] = {
+ 5.0121582440e-09, /* atan(0.5)lo 0x31ac3769 */
+ 3.7748947079e-08, /* atan(1.0)lo 0x33222168 */
+ 3.4473217170e-08, /* atan(1.5)lo 0x33140fb4 */
+ 7.5497894159e-08, /* atan(inf)lo 0x33a22168 */
+};
+
+static const float aT[] = {
+ 3.3333334327e-01, /* 0x3eaaaaaa */
+ -2.0000000298e-01, /* 0xbe4ccccd */
+ 1.4285714924e-01, /* 0x3e124925 */
+ -1.1111110449e-01, /* 0xbde38e38 */
+ 9.0908870101e-02, /* 0x3dba2e6e */
+ -7.6918758452e-02, /* 0xbd9d8795 */
+ 6.6610731184e-02, /* 0x3d886b35 */
+ -5.8335702866e-02, /* 0xbd6ef16b */
+ 4.9768779427e-02, /* 0x3d4bda59 */
+ -3.6531571299e-02, /* 0xbd15a221 */
+ 1.6285819933e-02, /* 0x3c8569d7 */
+};
+
+static const float
+huge = 1.0e+30,
+tiny = 1.0e-30,
+one = 1.0f;
+
+float atan_wrapper(float x)
{
- double w,s1,s2,z;
- int ix,hx,id;
+ float w,s1,s2,z;
+ int32_t ix,hx,id;
- hx = __HI(x);
+ GET_FLOAT_WORD(hx,x);
ix = hx&0x7fffffff;
- if(ix>=0x44100000) { /* if |x| >= 2^66 */
- if(ix>0x7ff00000||
- (ix==0x7ff00000&&(__LO(x)!=0)))
+ if(ix>=0x50800000) { /* if |x| >= 2^34 */
+ if(ix>0x7f800000)
return x+x; /* NaN */
if(hx>0) return atanhi[3]+atanlo[3];
else return -atanhi[3]-atanlo[3];
- } if (ix < 0x3fdc0000) { /* |x| < 0.4375 */
- if (ix < 0x3e200000) { /* |x| < 2^-29 */
+ } if (ix < 0x3ee00000) { /* |x| < 0.4375 */
+ if (ix < 0x31000000) { /* |x| < 2^-29 */
if(huge+x>one) return x; /* raise inexact */
}
id = -1;
} else {
- x = fabs(x);
- if (ix < 0x3ff30000) { /* |x| < 1.1875 */
- if (ix < 0x3fe60000) { /* 7/16 <=|x|<11/16 */
- id = 0; x = (2.0*x-one)/(2.0+x);
+ x = fabs_wrapper(x);
+ if (ix < 0x3f980000) { /* |x| < 1.1875 */
+ if (ix < 0x3f300000) { /* 7/16 <=|x|<11/16 */
+ id = 0; x = ((float)2.0*x-one)/((float)2.0+x);
} else { /* 11/16<=|x|< 19/16 */
id = 1; x = (x-one)/(x+one);
}
} else {
- if (ix < 0x40038000) { /* |x| < 2.4375 */
- id = 2; x = (x-1.5)/(one+1.5*x);
+ if (ix < 0x401c0000) { /* |x| < 2.4375 */
+ id = 2; x = (x-(float)1.5)/(one+(float)1.5*x);
} else { /* 2.4375 <= |x| < 2^66 */
- id = 3; x = -1.0/x;
+ id = 3; x = -(float)1.0/x;
}
}}
/* end of argument reduction */
@@ -312,444 +871,148 @@ double atan_wrapper(double x)
return (hx<0)? -z:z;
}
}
-/* @(#)e_sqrt.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.
- * ====================================================
- */
-/* __ieee754_sqrt(x)
- * Return correctly rounded sqrt.
- * ------------------------------------------
- * | Use the hardware sqrt if you have one |
- * ------------------------------------------
- * Method:
- * Bit by bit method using integer arithmetic. (Slow, but portable)
- * 1. Normalization
- * Scale x to y in [1,4) with even powers of 2:
- * find an integer k such that 1 <= (y=x*2^(2k)) < 4, then
- * sqrt(x) = 2^k * sqrt(y)
- * 2. Bit by bit computation
- * Let q = sqrt(y) truncated to i bit after binary point (q = 1),
- * i 0
- * i+1 2
- * s = 2*q , and y = 2 * ( y - q ). (1)
- * i i i i
- *
- * To compute q from q , one checks whether
- * i+1 i
- *
- * -(i+1) 2
- * (q + 2 ) <= y. (2)
- * i
- * -(i+1)
- * If (2) is false, then q = q ; otherwise q = q + 2 .
- * i+1 i i+1 i
- *
- * With some algebric manipulation, it is not difficult to see
- * that (2) is equivalent to
- * -(i+1)
- * s + 2 <= y (3)
- * i i
- *
- * The advantage of (3) is that s and y can be computed by
- * i i
- * the following recurrence formula:
- * if (3) is false
- *
- * s = s , y = y ; (4)
- * i+1 i i+1 i
- *
- * otherwise,
- * -i -(i+1)
- * s = s + 2 , y = y - s - 2 (5)
- * i+1 i i+1 i i
- *
- * One may easily use induction to prove (4) and (5).
- * Note. Since the left hand side of (3) contain only i+2 bits,
- * it does not necessary to do a full (53-bit) comparison
- * in (3).
- * 3. Final rounding
- * After generating the 53 bits result, we compute one more bit.
- * Together with the remainder, we can decide whether the
- * result is exact, bigger than 1/2ulp, or less than 1/2ulp
- * (it will never equal to 1/2ulp).
- * The rounding mode can be detected by checking whether
- * huge + tiny is equal to huge, and whether huge - tiny is
- * equal to huge for some floating point number "huge" and "tiny".
- *
- * Special cases:
- * sqrt(+-0) = +-0 ... exact
- * sqrt(inf) = inf
- * sqrt(-ve) = NaN ... with invalid signal
- * sqrt(NaN) = NaN ... with invalid signal for signaling NaN
- *
- * Other methods : see the appended file at the end of the program below.
- *---------------
- */
+/* Arc tangent from two variables, original. */
+
+static const float
+pi_o_4 = 7.8539818525e-01, /* 0x3f490fdb */
+pi_o_2 = 1.5707963705e+00, /* 0x3fc90fdb */
+pi = 3.1415927410e+00, /* 0x40490fdb */
+pi_lo = -8.7422776573e-08; /* 0xb3bbbd2e */
-double sqrt_wrapper(double x)
+float atan2_wrapper(float y, float x)
{
- double z;
- int sign = (int)0x80000000;
- unsigned r,t1,s1,ix1,q1;
- int ix0,s0,q,m,t,i;
+ float z;
+ int32_t k,m,hx,hy,ix,iy;
- ix0 = __HI(x); /* high word of x */
- ix1 = __LO(x); /* low word of x */
+ GET_FLOAT_WORD(hx,x);
+ ix = hx&0x7fffffff;
+ GET_FLOAT_WORD(hy,y);
+ iy = hy&0x7fffffff;
+ if((ix>0x7f800000)||
+ (iy>0x7f800000)) /* x or y is NaN */
+ return x+y;
+ if(hx==0x3f800000) return atan_wrapper(y); /* x=1.0 */
+ m = ((hy>>31)&1)|((hx>>30)&2); /* 2*sign(x)+sign(y) */
+
+ /* when y = 0 */
+ if(iy==0) {
+ switch(m) {
+ case 0:
+ case 1: return y; /* atan(+-0,+anything)=+-0 */
+ case 2: return pi+tiny;/* atan(+0,-anything) = pi */
+ case 3: return -pi-tiny;/* atan(-0,-anything) =-pi */
+ }
+ }
+ /* when x = 0 */
+ if(ix==0) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
+
+ /* when x is INF */
+ if(ix==0x7f800000) {
+ if(iy==0x7f800000) {
+ switch(m) {
+ case 0: return pi_o_4+tiny;/* atan(+INF,+INF) */
+ case 1: return -pi_o_4-tiny;/* atan(-INF,+INF) */
+ case 2: return (float)3.0*pi_o_4+tiny;/*atan(+INF,-INF)*/
+ case 3: return (float)-3.0*pi_o_4-tiny;/*atan(-INF,-INF)*/
+ }
+ } else {
+ switch(m) {
+ case 0: return zero ; /* atan(+...,+INF) */
+ case 1: return -zero ; /* atan(-...,+INF) */
+ case 2: return pi+tiny ; /* atan(+...,-INF) */
+ case 3: return -pi-tiny ; /* atan(-...,-INF) */
+ }
+ }
+ }
+ /* when y is INF */
+ if(iy==0x7f800000) return (hy<0)? -pi_o_2-tiny: pi_o_2+tiny;
+
+ /* compute y/x */
+ k = (iy-ix)>>23;
+ if(k > 60) z=pi_o_2+(float)0.5*pi_lo; /* |y/x| > 2**60 */
+ else if(hx<0&&k<-60) z=0.0; /* |y|/x < -2**60 */
+ else z=atan_wrapper(fabs_wrapper(y/x)); /* safe to do y/x */
+ switch (m) {
+ case 0: return z ; /* atan(+,+) */
+ case 1: {
+ uint32_t zh;
+ GET_FLOAT_WORD(zh,z);
+ SET_FLOAT_WORD(z,zh ^ 0x80000000);
+ }
+ return z ; /* atan(-,+) */
+ case 2: return pi-(z-pi_lo);/* atan(+,-) */
+ default: /* case 3 */
+ return (z-pi_lo)-pi;/* atan(-,-) */
+ }
+}
+
+/* Square root function, original. */
+float sqrt_wrapper(float x)
+{
+ float z;
+ int32_t sign = (int)0x80000000;
+ int32_t ix,s,q,m,t,i;
+ uint32_t r;
+
+ GET_FLOAT_WORD(ix,x);
/* take care of Inf and NaN */
- if((ix0&0x7ff00000)==0x7ff00000) {
+ if((ix&0x7f800000)==0x7f800000) {
return x*x+x; /* sqrt(NaN)=NaN, sqrt(+inf)=+inf
sqrt(-inf)=sNaN */
}
/* take care of zero */
- if(ix0<=0) {
- if(((ix0&(~sign))|ix1)==0) return x;/* sqrt(+-0) = +-0 */
- else if(ix0<0)
+ if(ix<=0) {
+ if((ix&(~sign))==0) return x;/* sqrt(+-0) = +-0 */
+ else if(ix<0)
return (x-x)/(x-x); /* sqrt(-ve) = sNaN */
}
/* normalize x */
- m = (ix0>>20);
+ m = (ix>>23);
if(m==0) { /* subnormal x */
- while(ix0==0) {
- m -= 21;
- ix0 |= (ix1>>11); ix1 <<= 21;
- }
- for(i=0;(ix0&0x00100000)==0;i++) ix0<<=1;
+ for(i=0;(ix&0x00800000)==0;i++) ix<<=1;
m -= i-1;
- ix0 |= (ix1>>(32-i));
- ix1 <<= i;
- }
- m -= 1023; /* unbias exponent */
- ix0 = (ix0&0x000fffff)|0x00100000;
- if(m&1){ /* odd m, double x to make it even */
- ix0 += ix0 + ((ix1&sign)>>31);
- ix1 += ix1;
}
+ m -= 127; /* unbias exponent */
+ ix = (ix&0x007fffff)|0x00800000;
+ if(m&1) /* odd m, double x to make it even */
+ ix += ix;
m >>= 1; /* m = [m/2] */
/* generate sqrt(x) bit by bit */
- ix0 += ix0 + ((ix1&sign)>>31);
- ix1 += ix1;
- q = q1 = s0 = s1 = 0; /* [q,q1] = sqrt(x) */
- r = 0x00200000; /* r = moving bit from right to left */
+ ix += ix;
+ q = s = 0; /* q = sqrt(x) */
+ r = 0x01000000; /* r = moving bit from right to left */
while(r!=0) {
- t = s0+r;
- if(t<=ix0) {
- s0 = t+r;
- ix0 -= t;
+ t = s+r;
+ if(t<=ix) {
+ s = t+r;
+ ix -= t;
q += r;
}
- ix0 += ix0 + ((ix1&sign)>>31);
- ix1 += ix1;
- r>>=1;
- }
-
- r = sign;
- while(r!=0) {
- t1 = s1+r;
- t = s0;
- if((t<ix0)||((t==ix0)&&(t1<=ix1))) {
- s1 = t1+r;
- if(((t1&sign)==sign)&&(s1&sign)==0) s0 += 1;
- ix0 -= t;
- if (ix1 < t1) ix0 -= 1;
- ix1 -= t1;
- q1 += r;
- }
- ix0 += ix0 + ((ix1&sign)>>31);
- ix1 += ix1;
+ ix += ix;
r>>=1;
}
/* use floating add to find out rounding direction */
- if((ix0|ix1)!=0) {
+ if(ix!=0) {
z = one-tiny; /* trigger inexact flag */
if (z>=one) {
z = one+tiny;
- if (q1==(unsigned)0xffffffff) { q1=0; q += 1;}
- else if (z>one) {
- if (q1==(unsigned)0xfffffffe) q+=1;
- q1+=2;
- } else
- q1 += (q1&1);
+ if (z>one)
+ q += 2;
+ else
+ q += (q&1);
}
}
- ix0 = (q>>1)+0x3fe00000;
- ix1 = q1>>1;
- if ((q&1)==1) ix1 |= sign;
- ix0 += (m <<20);
- __HI(z) = ix0;
- __LO(z) = ix1;
+ ix = (q>>1)+0x3f000000;
+ ix += (m <<23);
+ SET_FLOAT_WORD(z,ix);
return z;
}
-/*
-Other methods (use floating-point arithmetic)
--------------
-(This is a copy of a drafted paper by Prof W. Kahan
-and K.C. Ng, written in May, 1986)
-
- Two algorithms are given here to implement sqrt(x)
- (IEEE double precision arithmetic) in software.
- Both supply sqrt(x) correctly rounded. The first algorithm (in
- Section A) uses newton iterations and involves four divisions.
- The second one uses reciproot iterations to avoid division, but
- requires more multiplications. Both algorithms need the ability
- to chop results of arithmetic operations instead of round them,
- and the INEXACT flag to indicate when an arithmetic operation
- is executed exactly with no roundoff error, all part of the
- standard (IEEE 754-1985). The ability to perform shift, add,
- subtract and logical AND operations upon 32-bit words is needed
- too, though not part of the standard.
-
-A. sqrt(x) by Newton Iteration
-
- (1) Initial approximation
-
- Let x0 and x1 be the leading and the trailing 32-bit words of
- a floating point number x (in IEEE double format) respectively
-
- 1 11 52 ...widths
- ------------------------------------------------------
- x: |s| e | f |
- ------------------------------------------------------
- msb lsb msb lsb ...order
-
-
- ------------------------ ------------------------
- x0: |s| e | f1 | x1: | f2 |
- ------------------------ ------------------------
-
- By performing shifts and subtracts on x0 and x1 (both regarded
- as integers), we obtain an 8-bit approximation of sqrt(x) as
- follows.
-
- k := (x0>>1) + 0x1ff80000;
- y0 := k - T1[31&(k>>15)]. ... y ~ sqrt(x) to 8 bits
- Here k is a 32-bit integer and T1[] is an integer array containing
- correction terms. Now magically the floating value of y (y's
- leading 32-bit word is y0, the value of its trailing word is 0)
- approximates sqrt(x) to almost 8-bit.
-
- Value of T1:
- static int T1[32]= {
- 0, 1024, 3062, 5746, 9193, 13348, 18162, 23592,
- 29598, 36145, 43202, 50740, 58733, 67158, 75992, 85215,
- 83599, 71378, 60428, 50647, 41945, 34246, 27478, 21581,
- 16499, 12183, 8588, 5674, 3403, 1742, 661, 130,};
-
- (2) Iterative refinement
-
- Apply Heron's rule three times to y, we have y approximates
- sqrt(x) to within 1 ulp (Unit in the Last Place):
-
- y := (y+x/y)/2 ... almost 17 sig. bits
- y := (y+x/y)/2 ... almost 35 sig. bits
- y := y-(y-x/y)/2 ... within 1 ulp
-
-
- Remark 1.
- Another way to improve y to within 1 ulp is:
-
- y := (y+x/y) ... almost 17 sig. bits to 2*sqrt(x)
- y := y - 0x00100006 ... almost 18 sig. bits to sqrt(x)
-
- 2
- (x-y )*y
- y := y + 2* ---------- ...within 1 ulp
- 2
- 3y + x
-
-
- This formula has one division fewer than the one above; however,
- it requires more multiplications and additions. Also x must be
- scaled in advance to avoid spurious overflow in evaluating the
- expression 3y*y+x. Hence it is not recommended uless division
- is slow. If division is very slow, then one should use the
- reciproot algorithm given in section B.
-
- (3) Final adjustment
-
- By twiddling y's last bit it is possible to force y to be
- correctly rounded according to the prevailing rounding mode
- as follows. Let r and i be copies of the rounding mode and
- inexact flag before entering the square root program. Also we
- use the expression y+-ulp for the next representable floating
- numbers (up and down) of y. Note that y+-ulp = either fixed
- point y+-1, or multiply y by nextafter(1,+-inf) in chopped
- mode.
-
- I := FALSE; ... reset INEXACT flag I
- R := RZ; ... set rounding mode to round-toward-zero
- z := x/y; ... chopped quotient, possibly inexact
- If(not I) then { ... if the quotient is exact
- if(z=y) {
- I := i; ... restore inexact flag
- R := r; ... restore rounded mode
- return sqrt(x):=y.
- } else {
- z := z - ulp; ... special rounding
- }
- }
- i := TRUE; ... sqrt(x) is inexact
- If (r=RN) then z=z+ulp ... rounded-to-nearest
- If (r=RP) then { ... round-toward-+inf
- y = y+ulp; z=z+ulp;
- }
- y := y+z; ... chopped sum
- y0:=y0-0x00100000; ... y := y/2 is correctly rounded.
- I := i; ... restore inexact flag
- R := r; ... restore rounded mode
- return sqrt(x):=y.
-
- (4) Special cases
-
- Square root of +inf, +-0, or NaN is itself;
- Square root of a negative number is NaN with invalid signal.
-
-
-B. sqrt(x) by Reciproot Iteration
-
- (1) Initial approximation
-
- Let x0 and x1 be the leading and the trailing 32-bit words of
- a floating point number x (in IEEE double format) respectively
- (see section A). By performing shifs and subtracts on x0 and y0,
- we obtain a 7.8-bit approximation of 1/sqrt(x) as follows.
-
- k := 0x5fe80000 - (x0>>1);
- y0:= k - T2[63&(k>>14)]. ... y ~ 1/sqrt(x) to 7.8 bits
-
- Here k is a 32-bit integer and T2[] is an integer array
- containing correction terms. Now magically the floating
- value of y (y's leading 32-bit word is y0, the value of
- its trailing word y1 is set to zero) approximates 1/sqrt(x)
- to almost 7.8-bit.
-
- Value of T2:
- static int T2[64]= {
- 0x1500, 0x2ef8, 0x4d67, 0x6b02, 0x87be, 0xa395, 0xbe7a, 0xd866,
- 0xf14a, 0x1091b,0x11fcd,0x13552,0x14999,0x15c98,0x16e34,0x17e5f,
- 0x18d03,0x19a01,0x1a545,0x1ae8a,0x1b5c4,0x1bb01,0x1bfde,0x1c28d,
- 0x1c2de,0x1c0db,0x1ba73,0x1b11c,0x1a4b5,0x1953d,0x18266,0x16be0,
- 0x1683e,0x179d8,0x18a4d,0x19992,0x1a789,0x1b445,0x1bf61,0x1c989,
- 0x1d16d,0x1d77b,0x1dddf,0x1e2ad,0x1e5bf,0x1e6e8,0x1e654,0x1e3cd,
- 0x1df2a,0x1d635,0x1cb16,0x1be2c,0x1ae4e,0x19bde,0x1868e,0x16e2e,
- 0x1527f,0x1334a,0x11051,0xe951, 0xbe01, 0x8e0d, 0x5924, 0x1edd,};
-
- (2) Iterative refinement
-
- Apply Reciproot iteration three times to y and multiply the
- result by x to get an approximation z that matches sqrt(x)
- to about 1 ulp. To be exact, we will have
- -1ulp < sqrt(x)-z<1.0625ulp.
-
- ... set rounding mode to Round-to-nearest
- y := y*(1.5-0.5*x*y*y) ... almost 15 sig. bits to 1/sqrt(x)
- y := y*((1.5-2^-30)+0.5*x*y*y)... about 29 sig. bits to 1/sqrt(x)
- ... special arrangement for better accuracy
- z := x*y ... 29 bits to sqrt(x), with z*y<1
- z := z + 0.5*z*(1-z*y) ... about 1 ulp to sqrt(x)
-
- Remark 2. The constant 1.5-2^-30 is chosen to bias the error so that
- (a) the term z*y in the final iteration is always less than 1;
- (b) the error in the final result is biased upward so that
- -1 ulp < sqrt(x) - z < 1.0625 ulp
- instead of |sqrt(x)-z|<1.03125ulp.
-
- (3) Final adjustment
-
- By twiddling y's last bit it is possible to force y to be
- correctly rounded according to the prevailing rounding mode
- as follows. Let r and i be copies of the rounding mode and
- inexact flag before entering the square root program. Also we
- use the expression y+-ulp for the next representable floating
- numbers (up and down) of y. Note that y+-ulp = either fixed
- point y+-1, or multiply y by nextafter(1,+-inf) in chopped
- mode.
-
- R := RZ; ... set rounding mode to round-toward-zero
- switch(r) {
- case RN: ... round-to-nearest
- if(x<= z*(z-ulp)...chopped) z = z - ulp; else
- if(x<= z*(z+ulp)...chopped) z = z; else z = z+ulp;
- break;
- case RZ:case RM: ... round-to-zero or round-to--inf
- R:=RP; ... reset rounding mod to round-to-+inf
- if(x<z*z ... rounded up) z = z - ulp; else
- if(x>=(z+ulp)*(z+ulp) ...rounded up) z = z+ulp;
- break;
- case RP: ... round-to-+inf
- if(x>(z+ulp)*(z+ulp)...chopped) z = z+2*ulp; else
- if(x>z*z ...chopped) z = z+ulp;
- break;
- }
-
- Remark 3. The above comparisons can be done in fixed point. For
- example, to compare x and w=z*z chopped, it suffices to compare
- x1 and w1 (the trailing parts of x and w), regarding them as
- two's complement integers.
-
- ...Is z an exact square root?
- To determine whether z is an exact square root of x, let z1 be the
- trailing part of z, and also let x0 and x1 be the leading and
- trailing parts of x.
-
- If ((z1&0x03ffffff)!=0) ... not exact if trailing 26 bits of z!=0
- I := 1; ... Raise Inexact flag: z is not exact
- else {
- j := 1 - [(x0>>20)&1] ... j = logb(x) mod 2
- k := z1 >> 26; ... get z's 25-th and 26-th
- fraction bits
- I := i or (k&j) or ((k&(j+j+1))!=(x1&3));
- }
- R:= r ... restore rounded mode
- return sqrt(x):=z.
-
- If multiplication is cheaper then the foregoing red tape, the
- Inexact flag can be evaluated by
-
- I := i;
- I := (z*z!=x) or I.
-
- Note that z*z can overwrite I; this value must be sensed if it is
- True.
-
- Remark 4. If z*z = x exactly, then bit 25 to bit 0 of z1 must be
- zero.
-
- --------------------
- z1: | f2 |
- --------------------
- bit 31 bit 0
-
- Further more, bit 27 and 26 of z1, bit 0 and 1 of x1, and the odd
- or even of logb(x) have the following relations:
-
- -------------------------------------------------
- bit 27,26 of z1 bit 1,0 of x1 logb(x)
- -------------------------------------------------
- 00 00 odd and even
- 01 01 even
- 10 10 odd
- 10 00 even
- 11 01 even
- -------------------------------------------------
-
- (4) Special cases (see (4) of Section A).
-
- */
-
-
/* @(#)e_acos.c 1.3 95/01/18 */
/*
* ====================================================
@@ -1029,300 +1292,50 @@ uint64_t strtouq_wrapper(CONST char *nptr, char **endptr, int base)
return strtol(nptr, endptr, base);
}
-/*
- * ====================================================
- * Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
- *
- * Permission to use, copy, modify, and distribute this
- * software is freely granted, provided that this notice
- * is preserved.
- * ====================================================
- */
-
-/* __ieee754_pow(x,y) return x**y
- *
- * n
- * Method: Let x = 2 * (1+f)
- * 1. Compute and return log2(x) in two pieces:
- * log2(x) = w1 + w2,
- * where w1 has 53-24 = 29 bit trailing zeros.
- * 2. Perform y*log2(x) = n+y' by simulating muti-precision
- * arithmetic, where |y'|<=0.5.
- * 3. Return x**y = 2**n*exp(y'*log2)
- *
- * Special cases:
- * 1. (anything) ** 0 is 1
- * 2. (anything) ** 1 is itself
- * 3. (anything) ** NAN is NAN
- * 4. NAN ** (anything except 0) is NAN
- * 5. +-(|x| > 1) ** +INF is +INF
- * 6. +-(|x| > 1) ** -INF is +0
- * 7. +-(|x| < 1) ** +INF is +0
- * 8. +-(|x| < 1) ** -INF is +INF
- * 9. +-1 ** +-INF is NAN
- * 10. +0 ** (+anything except 0, NAN) is +0
- * 11. -0 ** (+anything except 0, NAN, odd integer) is +0
- * 12. +0 ** (-anything except 0, NAN) is +INF
- * 13. -0 ** (-anything except 0, NAN, odd integer) is +INF
- * 14. -0 ** (odd integer) = -( +0 ** (odd integer) )
- * 15. +INF ** (+anything except 0,NAN) is +INF
- * 16. +INF ** (-anything except 0,NAN) is +0
- * 17. -INF ** (anything) = -0 ** (-anything)
- * 18. (-anything) ** (integer) is (-1)**(integer)*(+anything**integer)
- * 19. (-anything except 0 and inf) ** (non-integer) is NAN
- *
- * Accuracy:
- * pow(x,y) returns x**y nearly rounded. In particular
- * pow(integer,integer)
- * always returns the correct integer provided it is
- * representable.
- *
- * Constants :
- * The hexadecimal values are the intended ones for the following
- * constants. The decimal values may be used, provided that the
- * compiler will convert from decimal to binary accurately enough
- * to produce the hexadecimal values shown.
- */
-
-double scalbn_wrapper (double x, int n);
-
-#ifdef __STDC__
-static const double
-#else
-static double
-#endif
-bp[] = {1.0, 1.5,},
-dp_h[] = { 0.0, 5.84962487220764160156e-01,}, /* 0x3FE2B803, 0x40000000 */
-dp_l[] = { 0.0, 1.35003920212974897128e-08,}, /* 0x3E4CFDEB, 0x43CFD006 */
-two = 2.0,
-two53 = 9007199254740992.0, /* 0x43400000, 0x00000000 */
- /* poly coefs for (3/2)*(log(x)-2s-2/3*s**3 */
-L1 = 5.99999999999994648725e-01, /* 0x3FE33333, 0x33333303 */
-L2 = 4.28571428578550184252e-01, /* 0x3FDB6DB6, 0xDB6FABFF */
-L3 = 3.33333329818377432918e-01, /* 0x3FD55555, 0x518F264D */
-L4 = 2.72728123808534006489e-01, /* 0x3FD17460, 0xA91D4101 */
-L5 = 2.30660745775561754067e-01, /* 0x3FCD864A, 0x93C9DB65 */
-L6 = 2.06975017800338417784e-01, /* 0x3FCA7E28, 0x4A454EEF */
-P1 = 1.66666666666666019037e-01, /* 0x3FC55555, 0x5555553E */
-P2 = -2.77777777770155933842e-03, /* 0xBF66C16C, 0x16BEBD93 */
-P3 = 6.61375632143793436117e-05, /* 0x3F11566A, 0xAF25DE2C */
-P4 = -1.65339022054652515390e-06, /* 0xBEBBBD41, 0xC5D26BF1 */
-P5 = 4.13813679705723846039e-08, /* 0x3E663769, 0x72BEA4D0 */
-lg2 = 6.93147180559945286227e-01, /* 0x3FE62E42, 0xFEFA39EF */
-lg2_h = 6.93147182464599609375e-01, /* 0x3FE62E43, 0x00000000 */
-lg2_l = -1.90465429995776804525e-09, /* 0xBE205C61, 0x0CA86C39 */
-ovt = 8.0085662595372944372e-0017, /* -(1024-log2(ovfl+.5ulp)) */
-cp = 9.61796693925975554329e-01, /* 0x3FEEC709, 0xDC3A03FD =2/(3ln2) */
-cp_h = 9.61796700954437255859e-01, /* 0x3FEEC709, 0xE0000000 =(float)cp */
-cp_l = -7.02846165095275826516e-09, /* 0xBE3E2FE0, 0x145B01F5 =tail of cp_h*/
-ivln2 = 1.44269504088896338700e+00, /* 0x3FF71547, 0x652B82FE =1/ln2 */
-ivln2_h = 1.44269502162933349609e+00, /* 0x3FF71547, 0x60000000 =24b 1/ln2*/
-ivln2_l = 1.92596299112661746887e-08; /* 0x3E54AE0B, 0xF85DDF44 =1/ln2 tail*/
-
-double pow_wrapper(double x, double y)
+/* Power function, taken from glibc-2.8 and dietlibc-0.32 */
+float pow_wrapper(float x, float y)
{
- double z,ax,z_h,z_l,p_h,p_l;
- double y1,t1,t2,r,s,t,u,v,w;
- int i0,i1,i,j,k,yisint,n;
- int hx,hy,ix,iy;
- unsigned lx,ly;
-
- i0 = ((*(int*)&one)>>29)^1; i1=1-i0;
- hx = __HI(x); lx = __LO(x);
- hy = __HI(y); ly = __LO(y);
- ix = hx&0x7fffffff; iy = hy&0x7fffffff;
-
- /* y==zero: x**0 = 1 */
- if((iy|ly)==0) return one;
-
- /* +-NaN return x+y */
- if(ix > 0x7ff00000 || ((ix==0x7ff00000)&&(lx!=0)) ||
- iy > 0x7ff00000 || ((iy==0x7ff00000)&&(ly!=0)))
- return x+y;
-
- /* determine if y is an odd int when x < 0
- * yisint = 0 ... y is not an integer
- * yisint = 1 ... y is an odd int
- * yisint = 2 ... y is an even int
- */
- yisint = 0;
- if(hx<0) {
- if(iy>=0x43400000) yisint = 2; /* even integer y */
- else if(iy>=0x3ff00000) {
- k = (iy>>20)-0x3ff; /* exponent */
- if(k>20) {
- j = ly>>(52-k);
- if((j<<(52-k))==ly) yisint = 2-(j&1);
- } else if(ly==0) {
- j = iy>>(20-k);
- if((j<<(20-k))==iy) yisint = 2-(j&1);
- }
- }
- }
-
- /* special value of y */
- if(ly==0) {
- if (iy==0x7ff00000) { /* y is +-inf */
- if(((ix-0x3ff00000)|lx)==0)
- return y - y; /* inf**+-1 is NaN */
- else if (ix >= 0x3ff00000)/* (|x|>1)**+-inf = inf,0 */
- return (hy>=0)? y: zero;
- else /* (|x|<1)**-,+inf = inf,0 */
- return (hy<0)?-y: zero;
- }
- if(iy==0x3ff00000) { /* y is +-1 */
- if(hy<0) return one/x; else return x;
- }
- if(hy==0x40000000) return x*x; /* y is 2 */
- if(hy==0x3fe00000) { /* y is 0.5 */
- if(hx>=0) /* x >= +0 */
- return sqrt(x);
- }
- }
+ unsigned int e;
+ float result;
+
+ /* Special cases 0^x */
+ if(x == 0.0f)
+ {
+ if(y > 0.0f)
+ return 0.0f;
+ else if(y == 0.0f)
+ return 1.0f;
+ else
+ return 1.0f / x;
+ }
- ax = fabs(x);
- /* special value of x */
- if(lx==0) {
- if(ix==0x7ff00000||ix==0||ix==0x3ff00000){
- z = ax; /*x is +-0,+-inf,+-1*/
- if(hy<0) z = one/z; /* z = (1/|x|) */
- if(hx<0) {
- if(((ix-0x3ff00000)|yisint)==0) {
- z = (z-z)/(z-z); /* (-1)**non-int is NaN */
- } else if(yisint==1)
- z = -z; /* (x<0)**odd = -(|x|**odd) */
- }
- return z;
- }
+ /* Special case x^n where n is integer */
+ if(y == (int) (e = (int) y))
+ {
+ if((int) e < 0)
+ {
+ e = -e;
+ x = 1.0f / x;
}
- n = (hx>>31)+1;
+ result = 1.0f;
- /* (x<0)**(non-int) is NaN */
- if((n|yisint)==0) return (x-x)/(x-x);
+ while(1)
+ {
+ if(e & 1)
+ result *= x;
- s = one; /* s (sign of result -ve**odd) = -1 else = 1 */
- if((n|(yisint-1))==0) s = -one;/* (-ve)**(odd int) */
+ if((e >>= 1) == 0)
+ break;
- /* |y| is huge */
- if(iy>0x41e00000) { /* if |y| > 2**31 */
- if(iy>0x43f00000){ /* if |y| > 2**64, must o/uflow */
- if(ix<=0x3fefffff) return (hy<0)? huge*huge:tiny*tiny;
- if(ix>=0x3ff00000) return (hy>0)? huge*huge:tiny*tiny;
- }
- /* over/underflow if x is not close to one */
- if(ix<0x3fefffff) return (hy<0)? s*huge*huge:s*tiny*tiny;
- if(ix>0x3ff00000) return (hy>0)? s*huge*huge:s*tiny*tiny;
- /* now |1-x| is tiny <= 2**-20, suffice to compute
- log(x) by x-x^2/2+x^3/3-x^4/4 */
- t = ax-one; /* t has 20 trailing zeros */
- w = (t*t)*(0.5-t*(0.3333333333333333333333-t*0.25));
- u = ivln2_h*t; /* ivln2_h has 21 sig. bits */
- v = t*ivln2_l-w*ivln2;
- t1 = u+v;
- __LO(t1) = 0;
- t2 = v-(t1-u);
- } else {
- double ss,s2,s_h,s_l,t_h,t_l;
- n = 0;
- /* take care subnormal number */
- if(ix<0x00100000)
- {ax *= two53; n -= 53; ix = __HI(ax); }
- n += ((ix)>>20)-0x3ff;
- j = ix&0x000fffff;
- /* determine interval */
- ix = j|0x3ff00000; /* normalize ix */
- if(j<=0x3988E) k=0; /* |x|<sqrt(3/2) */
- else if(j<0xBB67A) k=1; /* |x|<sqrt(3) */
- else {k=0;n+=1;ix -= 0x00100000;}
- __HI(ax) = ix;
-
- /* compute ss = s_h+s_l = (x-1)/(x+1) or (x-1.5)/(x+1.5) */
- u = ax-bp[k]; /* bp[0]=1.0, bp[1]=1.5 */
- v = one/(ax+bp[k]);
- ss = u*v;
- s_h = ss;
- __LO(s_h) = 0;
- /* t_h=ax+bp[k] High */
- t_h = zero;
- __HI(t_h)=((ix>>1)|0x20000000)+0x00080000+(k<<18);
- t_l = ax - (t_h-bp[k]);
- s_l = v*((u-s_h*t_h)-s_h*t_l);
- /* compute log(ax) */
- s2 = ss*ss;
- r = s2*s2*(L1+s2*(L2+s2*(L3+s2*(L4+s2*(L5+s2*L6)))));
- r += s_l*(s_h+ss);
- s2 = s_h*s_h;
- t_h = 3.0+s2+r;
- __LO(t_h) = 0;
- t_l = r-((t_h-3.0)-s2);
- /* u+v = ss*(1+...) */
- u = s_h*t_h;
- v = s_l*t_h+t_l*ss;
- /* 2/(3log2)*(ss+...) */
- p_h = u+v;
- __LO(p_h) = 0;
- p_l = v-(p_h-u);
- z_h = cp_h*p_h; /* cp_h+cp_l = 2/(3*log2) */
- z_l = cp_l*p_h+p_l*cp+dp_l[k];
- /* log2(ax) = (ss+..)*2/(3*log2) = n + dp_h + z_h + z_l */
- t = (double)n;
- t1 = (((z_h+z_l)+dp_h[k])+t);
- __LO(t1) = 0;
- t2 = z_l-(((t1-t)-dp_h[k])-z_h);
+ x *= x;
}
- /* split up y into y1+y2 and compute (y1+y2)*(t1+t2) */
- y1 = y;
- __LO(y1) = 0;
- p_l = (y-y1)*t1+y*t2;
- p_h = y1*t1;
- z = p_l+p_h;
- j = __HI(z);
- i = __LO(z);
- if (j>=0x40900000) { /* z >= 1024 */
- if(((j-0x40900000)|i)!=0) /* if z > 1024 */
- return s*huge*huge; /* overflow */
- else {
- if(p_l+ovt>z-p_h) return s*huge*huge; /* overflow */
- }
- } else if((j&0x7fffffff)>=0x4090cc00 ) { /* z <= -1075 */
- if(((j-0xc090cc00)|i)!=0) /* z < -1075 */
- return s*tiny*tiny; /* underflow */
- else {
- if(p_l<=z-p_h) return s*tiny*tiny; /* underflow */
- }
- }
- /*
- * compute 2**(p_h+p_l)
- */
- i = j&0x7fffffff;
- k = (i>>20)-0x3ff;
- n = 0;
- if(i>0x3fe00000) { /* if |z| > 0.5, set n = [z+0.5] */
- n = j+(0x00100000>>(k+1));
- k = ((n&0x7fffffff)>>20)-0x3ff; /* new k for n */
- t = zero;
- __HI(t) = (n&~(0x000fffff>>k));
- n = ((n&0x000fffff)|0x00100000)>>(20-k);
- if(j<0) n = -n;
- p_h -= t;
- }
- t = p_l+p_h;
- __LO(t) = 0;
- u = t*lg2_h;
- v = (p_l-(t-p_h))*lg2+t*lg2_l;
- z = u+v;
- w = v-(z-u);
- t = z*z;
- t1 = z - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))));
- r = (z*t1)/(t1-two)-(w+z*w);
- z = one-(r-z);
- j = __HI(z);
- j += (n<<20);
- if((j>>20)<=0) z = scalbn_wrapper(z,n); /* subnormal output */
- else __HI(z) += (n<<20);
- return s*z;
+ return result;
+ }
+
+ /* Normal case */
+ return rb_exp(rb_log(x) * y);
}
/* @(#)s_copysign.c 1.3 95/01/18 */
@@ -1403,119 +1416,146 @@ double scalbn_wrapper (double x, int n)
return x*twom54;
}
-/* @(#)s_ceil.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.
- * ====================================================
- */
+/* horrible hack */
+float ceil_wrapper(float x)
+{
+ return floor_wrapper(x) + 1.0;
+}
-/*
- * ceil(x)
- * Return x rounded toward -inf to integral value
- * Method:
- * Bit twiddling.
- * Exception:
- * Inexact flag raised if x not equal to ceil(x).
- */
+/* Implementation of strtod() and atof(),
+ taken from SanOS (http://www.jbox.dk/sanos/). */
+static int rb_errno = 0;
-double ceil_wrapper(double x)
+static double rb_strtod(const char *str, char **endptr)
{
- int i0,i1,j0;
- unsigned i,j;
- i0 = __HI(x);
- i1 = __LO(x);
- j0 = ((i0>>20)&0x7ff)-0x3ff;
- if(j0<20) {
- if(j0<0) { /* raise inexact if x != 0 */
- if(huge+x>0.0) {/* return 0*sign(x) if |x|<1 */
- if(i0<0) {i0=0x80000000;i1=0;}
- else if((i0|i1)!=0) { i0=0x3ff00000;i1=0;}
- }
- } else {
- i = (0x000fffff)>>j0;
- if(((i0&i)|i1)==0) return x; /* x is integral */
- if(huge+x>0.0) { /* raise inexact flag */
- if(i0>0) i0 += (0x00100000)>>j0;
- i0 &= (~i); i1=0;
- }
- }
- } else if (j0>51) {
- if(j0==0x400) return x+x; /* inf or NaN */
- else return x; /* x is integral */
- } else {
- i = ((unsigned)(0xffffffff))>>(j0-20);
- if((i1&i)==0) return x; /* x is integral */
- if(huge+x>0.0) { /* raise inexact flag */
- if(i0>0) {
- if(j0==20) i0+=1;
- else {
- j = i1 + (1<<(52-j0));
- if(j<i1) i0+=1; /* got a carry */
- i1 = j;
- }
- }
- i1 &= (~i);
- }
- }
- __HI(x) = i0;
- __LO(x) = i1;
- return x;
-}
+ double number;
+ int exponent;
+ int negative;
+ char *p = (char *) str;
+ double p10;
+ int n;
+ int num_digits;
+ int num_decimals;
+
+ /* Reset Rockbox errno -- W.B. */
+#ifdef ROCKBOX
+ rb_errno = 0;
+#endif
+
+ // Skip leading whitespace
+ while (isspace(*p)) p++;
+
+ // Handle optional sign
+ negative = 0;
+ switch (*p)
+ {
+ case '-': negative = 1; // Fall through to increment position
+ case '+': p++;
+ }
+
+ number = 0.;
+ exponent = 0;
+ num_digits = 0;
+ num_decimals = 0;
+
+ // Process string of digits
+ while (isdigit(*p))
+ {
+ number = number * 10. + (*p - '0');
+ p++;
+ num_digits++;
+ }
+
+ // Process decimal part
+ if (*p == '.')
+ {
+ p++;
+
+ while (isdigit(*p))
+ {
+ number = number * 10. + (*p - '0');
+ p++;
+ num_digits++;
+ num_decimals++;
+ }
+
+ exponent -= num_decimals;
+ }
-/* originally from python */
+ if (num_digits == 0)
+ {
+#ifdef ROCKBOX
+ rb_errno = 1;
+#else
+ errno = ERANGE;
+#endif
+ return 0.0;
+ }
+
+ // Correct for sign
+ if (negative) number = -number;
+
+ // Process an exponent string
+ if (*p == 'e' || *p == 'E')
+ {
+ // Handle optional sign
+ negative = 0;
+ switch(*++p)
+ {
+ case '-': negative = 1; // Fall through to increment pos
+ case '+': p++;
+ }
+
+ // Process string of digits
+ n = 0;
+ while (isdigit(*p))
+ {
+ n = n * 10 + (*p - '0');
+ p++;
+ }
+
+ if (negative)
+ exponent -= n;
+ else
+ exponent += n;
+ }
+
+#ifndef ROCKBOX
+ if (exponent < DBL_MIN_EXP || exponent > DBL_MAX_EXP)
+ {
+ errno = ERANGE;
+ return HUGE_VAL;
+ }
+#endif
+
+ // Scale the result
+ p10 = 10.;
+ n = exponent;
+ if (n < 0) n = -n;
+ while (n)
+ {
+ if (n & 1)
+ {
+ if (exponent < 0)
+ number /= p10;
+ else
+ number *= p10;
+ }
+ n >>= 1;
+ p10 *= p10;
+ }
-/* Just in case you haven't got an atof() around...
- This one doesn't check for bad syntax or overflow,
- and is slow and inaccurate.
- But it's good enough for the occasional string literal... */
+#ifndef ROCKBOX
+ if (number == HUGE_VAL) errno = ERANGE;
+#endif
+ if (endptr) *endptr = p;
-#include <ctype.h>
+ return number;
+}
-double atof_wrapper(char *s)
+double atof_wrapper(const char *str)
{
- double a = 0.0;
- int e = 0;
- int c;
- while ((c = *s++) != '\0' && isdigit(c)) {
- a = a*10.0 + (c - '0');
- }
- if (c == '.') {
- while ((c = *s++) != '\0' && isdigit(c)) {
- a = a*10.0 + (c - '0');
- e = e-1;
- }
- }
- if (c == 'e' || c == 'E') {
- int sign = 1;
- int i = 0;
- c = *s++;
- if (c == '+')
- c = *s++;
- else if (c == '-') {
- c = *s++;
- sign = -1;
- }
- while (isdigit(c)) {
- i = i*10 + (c - '0');
- c = *s++;
- }
- e += i*sign;
- }
- while (e > 0) {
- a *= 10.0;
- e--;
- }
- while (e < 0) {
- a *= 0.1;
- e++;
- }
- return a;
+ return rb_strtod(str, NULL);
}
/*
diff --git a/apps/plugins/puzzles/rockbox.c b/apps/plugins/puzzles/rockbox.c
index 32fc38e..c03088f 100644
--- a/apps/plugins/puzzles/rockbox.c
+++ b/apps/plugins/puzzles/rockbox.c
@@ -27,7 +27,12 @@
#include "lib/xlcd.h"
/* how many ticks between timer callbacks */
-#define TIMER_INTERVAL (HZ / 33)
+#define TIMER_INTERVAL (HZ / 50)
+#define BG_R .9f
+#define BG_G .9f
+#define BG_B .9f
+
+#define BG_COLOR LCD_RGBPACK((int)(255*BG_R), (int)(255*BG_R), (int)(255*BG_B))
static midend *me = NULL;
static unsigned *colors = NULL;
@@ -90,6 +95,7 @@ static void rb_color(int n)
static void rb_draw_text(void *handle, int x, int y, int fonttype,
int fontsize, int align, int colour, char *text)
{
+ (void) fontsize;
LOGF("rb_draw_text(%d %d %s)", x, y, text);
offset_coords(&x, &y);
@@ -148,13 +154,112 @@ static void rb_draw_line(void *handle, int x1, int y1, int x2, int y2,
rb->lcd_drawline(x1, y1, x2, y2);
}
+/*
+ * draw filled polygon
+ * originally by Sebastian Leonhardt (ulmutul)
+ * 'count' : number of coordinate pairs
+ * 'pxy': array of coordinates. pxy[0]=x0,pxy[1]=y0,...
+ * note: provide space for one extra coordinate, because the starting point
+ * will automatically be inserted as end point.
+ */
+
+/*
+ * helper function:
+ * find points of intersection between polygon and scanline
+ */
+
+#define MAX_INTERSECTION 32
+
+static void fill_poly_line(int scanline, int count, int *pxy)
+{
+ int i;
+ int j;
+ int num_of_intersects;
+ int direct, old_direct;
+ //intersections of every line with scanline (y-coord)
+ int intersection[MAX_INTERSECTION];
+ /* add starting point as ending point */
+ pxy[count*2] = pxy[0];
+ pxy[count*2+1] = pxy[1];
+
+ old_direct=0;
+ num_of_intersects=0;
+ for (i=0; i<count*2; i+=2) {
+ int x1=pxy[i];
+ int y1=pxy[i+1];
+ int x2=pxy[i+2];
+ int y2=pxy[i+3];
+ // skip if line is outside of scanline
+ if (y1 < y2) {
+ if (scanline < y1 || scanline > y2)
+ continue;
+ }
+ else {
+ if (scanline < y2 || scanline > y1)
+ continue;
+ }
+ // calculate x-coord of intersection
+ if (y1==y2) {
+ direct=0;
+ }
+ else {
+ direct = y1>y2 ? 1 : -1;
+ // omit double intersections, if both lines lead in the same direction
+ intersection[num_of_intersects] =
+ x1+((scanline-y1)*(x2-x1))/(y2-y1);
+ if ( (direct!=old_direct)
+ || (intersection[num_of_intersects] != intersection[num_of_intersects-1])
+ )
+ ++num_of_intersects;
+ }
+ old_direct = direct;
+ }
+
+ // sort points of intersection
+ for (i=0; i<num_of_intersects-1; ++i) {
+ for (j=i+1; j<num_of_intersects; ++j) {
+ if (intersection[j]<intersection[i]) {
+ int temp=intersection[i];
+ intersection[i]=intersection[j];
+ intersection[j]=temp;
+ }
+ }
+ }
+ // draw
+ for (i=0; i<num_of_intersects; i+=2) {
+ rb->lcd_hline(intersection[i], intersection[i+1], scanline);
+ }
+}
+
+/* two extra elements at end of pxy needed */
+static void v_fillarea(int count, int *pxy)
+{
+ int i;
+ int y1, y2;
+
+ // find min and max y coords
+ y1=y2=pxy[1];
+ for (i=3; i<count*2; i+=2) {
+ if (pxy[i] < y1) y1 = pxy[i];
+ else if (pxy[i] > y2) y2 = pxy[i];
+ }
+
+ for (i=y1; i<=y2; ++i) {
+ fill_poly_line(i, count, pxy);
+ }
+}
+
static void rb_draw_poly(void *handle, int *coords, int npoints,
int fillcolour, int outlinecolour)
{
LOGF("rb_draw_poly");
+
if(fillcolour >= 0)
{
rb_color(fillcolour);
+#if 1
+ /* serious hack: draw a bunch of triangles between adjacent points */
+ /* this generally works, even with some concave polygons */
for(int i = 2; i < npoints; ++i)
{
int x1, y1, x2, y2, x3, y3;
@@ -170,9 +275,36 @@ static void rb_draw_poly(void *handle, int *coords, int npoints,
xlcd_filltriangle(x1, y1,
x2, y2,
x3, y3);
+
+#if 0
+ rb->lcd_set_foreground(LCD_RGBPACK(255,0,0));
+ rb->lcd_drawpixel(x1, y1);
+ rb->lcd_drawpixel(x2, y2);
+ rb->lcd_drawpixel(x3, y3);
+ rb->lcd_update();
+ rb->sleep(HZ);
+ rb_color(fillcolour);
+ rb->lcd_drawpixel(x1, y1);
+ rb->lcd_drawpixel(x2, y2);
+ rb->lcd_drawpixel(x3, y3);
+ rb->lcd_update();
+#endif
+ }
+#else
+ int *pxy = smalloc(sizeof(int) * 2 * npoints + 2);
+ /* copy points, offsetted */
+ for(int i = 0; i < npoints; ++i)
+ {
+ pxy[2 * i + 0] = coords[2 * i + 0];
+ pxy[2 * i + 1] = coords[2 * i + 1];
+ offset_coords(&pxy[2*i+0], &pxy[2*i+1]);
}
+ v_fillarea(npoints, pxy);
+ sfree(pxy);
+#endif
}
+ /* draw outlines last so they're not covered by the fill */
assert(outlinecolour >= 0);
rb_color(outlinecolour);
@@ -187,6 +319,8 @@ static void rb_draw_poly(void *handle, int *coords, int npoints,
offset_coords(&x2, &y2);
rb->lcd_drawline(x1, y1,
x2, y2);
+ //rb->lcd_update();
+ //rb->sleep(HZ/2);
}
int x1, y1, x2, y2;
@@ -206,14 +340,16 @@ static void rb_draw_circle(void *handle, int cx, int cy, int radius,
{
LOGF("rb_draw_circle(%d, %d, %d)", cx, cy, radius);
offset_coords(&cx, &cy);
+
if(fillcolour >= 0)
{
rb_color(fillcolour);
- xlcd_fillcircle(cx, cy, radius);
+ xlcd_fillcircle(cx, cy, radius - 1);
}
+
assert(outlinecolour >= 0);
rb_color(outlinecolour);
- xlcd_drawcircle(cx, cy, radius);
+ xlcd_drawcircle(cx, cy, radius - 1);
}
struct blitter {
@@ -238,7 +374,6 @@ static void rb_blitter_free(void *handle, blitter *bl)
LOGF("rb_blitter_free");
sfree(bl->bmp.data);
sfree(bl);
- bl->bmp.data = NULL;
return;
}
@@ -270,6 +405,7 @@ static void trim_rect(int *x, int *y, int *w, int *h)
*h = y1 - y0;
}
+/* copy a section of the framebuffer */
static void rb_blitter_save(void *handle, blitter *bl, int x, int y)
{
/* no viewport offset */
@@ -284,7 +420,7 @@ static void rb_blitter_save(void *handle, blitter *bl, int x, int y)
for(int i = 0; i < h; ++i)
{
/* copy line-by-line */
- rb->memcpy(bl->bmp.data + sizeof(fb_data) * (y + i) * w,
+ rb->memcpy(bl->bmp.data + sizeof(fb_data) * i * w,
rb->lcd_framebuffer + (y + i) * LCD_WIDTH + x,
w * sizeof(fb_data));
}
@@ -307,7 +443,7 @@ static void rb_blitter_load(void *handle, blitter *bl, int x, int y)
trim_rect(&x, &y, &w, &h);
offset_coords(&x, &y);
- rb->lcd_bitmap(bl->bmp.data, x, y, w, h);
+ rb->lcd_bitmap((fb_data*)bl->bmp.data, x, y, w, h);
}
static void rb_draw_update(void *handle, int x, int y, int w, int h)
@@ -321,11 +457,44 @@ static void rb_end_draw(void *handle)
LOGF("rb_end_draw");
}
+static char *titlebar = NULL;
+
static void rb_status_bar(void *handle, char *text)
{
+ if(titlebar)
+ sfree(titlebar);
+ titlebar = dupstr(text);
LOGF("game title is %s\n", text);
}
+static void draw_title(void)
+{
+ const char *str = NULL;
+ if(titlebar)
+ str = titlebar;
+ else
+ str = midend_which_game(me)->name;
+
+ /* quick hack */
+ bool orig_clipped = clipped;
+ if(orig_clipped)
+ rb_unclip(NULL);
+
+ int h;
+ rb->lcd_setfont(FONT_UI);
+ rb->lcd_getstringsize(str, NULL, &h);
+
+ rb->lcd_set_foreground(BG_COLOR);
+ rb->lcd_fillrect(0, LCD_HEIGHT - h, LCD_WIDTH, h);
+
+ rb->lcd_set_foreground(LCD_BLACK);
+ rb->lcd_putsxy(0, LCD_HEIGHT - h, str);
+ rb->lcd_update_rect(0, LCD_HEIGHT - h, LCD_WIDTH, h);
+
+ if(orig_clipped)
+ rb_clip(NULL, clip_rect.x, clip_rect.y, clip_rect.width, clip_rect.height);
+}
+
static char *rb_text_fallback(void *handle, const char *const *strings,
int nstrings)
{
@@ -356,9 +525,9 @@ const drawing_api rb_drawing = {
void frontend_default_colour(frontend *fe, float *out)
{
- *out++ = 1.0;
- *out++ = 1.0;
- *out++ = 1.0;
+ *out++ = 0.9;
+ *out++ = 0.9;
+ *out++ = 0.9;
}
void fatal(char *fmt, ...)
@@ -381,6 +550,8 @@ void get_random_seed(void **randseed, int *randseedsize)
*randseed = snew(long);
long seed = *rb->current_tick;
rb->memcpy(*randseed, &seed, sizeof(seed));
+ //*(long*)*randseed = 42; // debug
+ //rb->splash(HZ, "DEBUG SEED ON");
*randseedsize = sizeof(long);
}
@@ -412,14 +583,14 @@ static int list_choose(const char *list_str, const char *title)
struct gui_synclist list;
- rb->gui_synclist_init(&list, &config_choices_formatter, list_str, false, 1, NULL);
+ rb->gui_synclist_init(&list, &config_choices_formatter, (void*)list_str, false, 1, NULL);
rb->gui_synclist_set_icon_callback(&list, NULL);
rb->gui_synclist_set_nb_items(&list, n);
rb->gui_synclist_limit_scroll(&list, false);
rb->gui_synclist_select_item(&list, 0);
- rb->gui_synclist_set_title(&list, title, NOICON);
+ rb->gui_synclist_set_title(&list, (char*)title, NOICON);
while (1)
{
rb->gui_synclist_draw(&list);
@@ -555,6 +726,16 @@ done:
free_cfg(config);
}
+static void quick_help(void)
+{
+ /* TODO */
+}
+
+static void full_help(void)
+{
+ /* TODO */
+}
+
static int pause_menu(void)
{
MENUITEM_STRINGLIST(menu, "Paused", NULL,
@@ -563,6 +744,8 @@ static int pause_menu(void)
"Undo",
"Redo",
"Solve",
+ "Quick Help",
+ "Extensive Help",
"Configure Game",
"Select Another Game",
"Quit");
@@ -588,7 +771,7 @@ static int pause_menu(void)
quit = true;
break;
case 3:
- if(!midend_can_undo(me))
+ if(!midend_can_redo(me))
rb->splash(HZ, "Cannot redo.");
else
midend_process_key(me, 0, 0, 'r');
@@ -603,17 +786,23 @@ static int pause_menu(void)
break;
}
case 5:
- config_menu();
+ quick_help();
break;
case 6:
- return -1;
+ full_help();
+ break;
case 7:
+ config_menu();
+ break;
+ case 8:
+ return -1;
+ case 9:
return -2;
default:
break;
}
}
- rb->lcd_set_background(LCD_WHITE);
+ rb->lcd_set_background(BG_COLOR);
rb->lcd_clear_display();
rb->lcd_update();
midend_force_redraw(me);
@@ -626,7 +815,17 @@ static int process_input(int tmo)
{
int state = 0;
static const struct button_mapping *plugin_contexts[] = { pla_main_ctx };
+
+#ifdef HAVE_ADJUSTABLE_CPU_FREQ
+ rb->cpu_boost(false); /* about to block for button input */
+#endif
+
int button = pluginlib_getaction(tmo, plugin_contexts, ARRAYLEN(plugin_contexts));
+
+#ifdef HAVE_ADJUSTABLE_CPU_FREQ
+ rb->cpu_boost(true);
+#endif
+
switch(button)
{
case PLA_UP:
@@ -645,8 +844,19 @@ static int process_input(int tmo)
state |= CURSOR_SELECT;
break;
case PLA_CANCEL:
+ {
want_redraw = false;
- return pause_menu();
+ /* quick hack to preserve the clipping state */
+ bool orig_clipped = clipped;
+ if(orig_clipped)
+ rb_unclip(NULL);
+
+ int rc = pause_menu();
+
+ if(orig_clipped)
+ rb_clip(NULL, clip_rect.x, clip_rect.y, clip_rect.width, clip_rect.height);
+ return rc;
+ }
default:
exit_on_usb(button);
break;
@@ -687,12 +897,20 @@ const char *formatter(char *buf, size_t n, int i, const char *unit)
static void fix_size(void)
{
- int w = LCD_WIDTH, h = LCD_HEIGHT;
+ int w = LCD_WIDTH, h = LCD_HEIGHT, h_x;
+ rb->lcd_setfont(FONT_UI);
+ rb->lcd_getstringsize("X", NULL, &h_x);
+ h -= h_x;
midend_size(me, &w, &h, TRUE);
}
static void init_for_game(const game *gm)
{
+ /* reset tlsf by nuking the signature */
+ /* will make any already-allocated memory point to garbage */
+ memset(giant_buffer, 0, 4);
+ init_memory_pool(sizeof(giant_buffer), giant_buffer);
+
me = midend_new(NULL, gm, &rb_drawing, NULL);
midend_new_game(me);
@@ -713,25 +931,38 @@ static void init_for_game(const game *gm)
}
/* seems to crash */
+ /* actually it doesn't matter if this memory is leaked as
+ * resetting tlsf makes it forget about it anyway */
//sfree(floatcolors);
rb->lcd_set_viewport(NULL);
rb->lcd_set_backdrop(NULL);
rb->lcd_set_foreground(LCD_BLACK);
- rb->lcd_set_background(LCD_WHITE);
+ rb->lcd_set_background(BG_COLOR);
rb->lcd_clear_display();
rb->lcd_update();
midend_force_redraw(me);
+ draw_title();
+}
+
+static void exit_handler(void)
+{
+#ifdef HAVE_ADJUSTABLE_CPU_FREQ
+ rb->cpu_boost(false);
+#endif
}
enum plugin_status plugin_start(const void *param)
{
(void) param;
- /* reset tlsf by nuking the signature */
- /* will make any already-allocated memory point to garbage */
- memset(giant_buffer, 0, 4);
- init_memory_pool(sizeof(giant_buffer), giant_buffer);
+ rb_atexit(exit_handler);
+#ifdef HAVE_ADJUSTABLE_CPU_FREQ
+ rb->cpu_boost(true);
+#endif
+
+ LOGF("acos(.5) = %f", acos(.5));
+ LOGF("sqrt(3)/2 = sin(60) = %f = %f", sqrt(3)/2, sin(PI/3));
int gm = 0;
while(1)
@@ -744,12 +975,21 @@ enum plugin_status plugin_start(const void *param)
while(1)
{
want_redraw = true;
+
+ draw_title();
+
int button = process_input(timer_on ? TIMER_INTERVAL : -1);
+
if(button < 0)
{
rb_unclip(NULL);
deactivate_timer(NULL);
midend_free(me);
+ if(titlebar)
+ {
+ sfree(titlebar);
+ titlebar = NULL;
+ }
/* new game */
if(button == -1)
break;
@@ -762,6 +1002,7 @@ enum plugin_status plugin_start(const void *param)
if(button)
midend_process_key(me, 0, 0, button);
+
if(want_redraw)
midend_redraw(me);
generated by cgit v1.1 at 2026-04-25 10:24:41 +0000