|
@@ -128,6 +128,8 @@
|
|
|
#include "math_libm.h"
|
|
|
#include "math_private.h"
|
|
|
|
|
|
+#include "SDL_assert.h"
|
|
|
+
|
|
|
static const int init_jk[] = {2,3,4,6}; /* initial value for jk */
|
|
|
|
|
|
static const double PIo2[] = {
|
|
@@ -147,13 +149,19 @@ one = 1.0,
|
|
|
two24 = 1.67772160000000000000e+07, /* 0x41700000, 0x00000000 */
|
|
|
twon24 = 5.96046447753906250000e-08; /* 0x3E700000, 0x00000000 */
|
|
|
|
|
|
-int32_t attribute_hidden __kernel_rem_pio2(double *x, double *y, int e0, int nx, int prec, const int32_t *ipio2)
|
|
|
+int32_t attribute_hidden __kernel_rem_pio2(double *x, double *y, int e0, int nx, const unsigned int prec, const int32_t *ipio2)
|
|
|
{
|
|
|
int32_t jz,jx,jv,jp,jk,carry,n,iq[20],i,j,k,m,q0,ih;
|
|
|
double z,fw,f[20],fq[20],q[20];
|
|
|
|
|
|
+ if (nx < 1) {
|
|
|
+ return 0;
|
|
|
+ }
|
|
|
+
|
|
|
/* initialize jk*/
|
|
|
+ SDL_assert(prec < SDL_arraysize(init_jk));
|
|
|
jk = init_jk[prec];
|
|
|
+ SDL_assert(jk > 0);
|
|
|
jp = jk;
|
|
|
|
|
|
/* determine jx,jv,q0, note that 3>q0 */
|
|
@@ -164,6 +172,7 @@ int32_t attribute_hidden __kernel_rem_pio2(double *x, double *y, int e0, int nx,
|
|
|
/* set up f[0] to f[jx+jk] where f[jx+jk] = ipio2[jv+jk] */
|
|
|
j = jv-jx; m = jx+jk;
|
|
|
for(i=0;i<=m;i++,j++) f[i] = (j<0)? zero : (double) ipio2[j];
|
|
|
+ for(i=m+1;i<SDL_arraysize(f);i++) f[i] = zero;
|
|
|
|
|
|
/* compute q[0],q[1],...q[jk] */
|
|
|
for (i=0;i<=jk;i++) {
|
|
@@ -179,6 +188,7 @@ recompute:
|
|
|
iq[i] = (int32_t)(z-two24*fw);
|
|
|
z = q[j-1]+fw;
|
|
|
}
|
|
|
+ iq[jz] = 0;
|
|
|
|
|
|
/* compute n */
|
|
|
z = scalbn(z,q0); /* actual value of z */
|
|
@@ -238,7 +248,8 @@ recompute:
|
|
|
/* chop off zero terms */
|
|
|
if(z==0.0) {
|
|
|
jz -= 1; q0 -= 24;
|
|
|
- while(iq[jz]==0) { jz--; q0-=24;}
|
|
|
+ SDL_assert(jz >= 0);
|
|
|
+ while(iq[jz]==0) { jz--; SDL_assert(jz >= 0); q0-=24;}
|
|
|
} else { /* break z into 24-bit if necessary */
|
|
|
z = scalbn(z,-q0);
|
|
|
if(z>=two24) {
|