SDL  2.0
k_rem_pio2.c File Reference
#include "math_libm.h"
#include "math_private.h"
#include "SDL_assert.h"
+ Include dependency graph for k_rem_pio2.c:

Go to the source code of this file.

Functions

 libm_hidden_proto (scalbn)
 
int attribute_hidden __kernel_rem_pio2 (x, y, int e0, int nx, int prec, ipio2)
 

Variables

static double PIo2 []
 
static double zero = 0.0
 
static double one = 1.0
 
static double two24 = 1.67772160000000000000e+07
 
static double twon24 = 5.96046447753906250000e-08
 

Function Documentation

◆ __kernel_rem_pio2()

int attribute_hidden __kernel_rem_pio2 ( x  ,
y  ,
int  e0,
int  nx,
int  prec,
ipio2   
)

Definition at line 176 of file k_rem_pio2.c.

References floor, i, j, k, one, scalbn, SDL_arraysize, SDL_assert, two24, twon24, and zero.

Referenced by __ieee754_rem_pio2().

181 {
182  int32_t jz, jx, jv, jp, jk, carry, n, iq[20], i, j, k, m, q0, ih;
183  double z, fw, f[20], fq[20], q[20];
184 
185  /* initialize jk */
186  SDL_assert((prec >= 0) && (prec < SDL_arraysize(init_jk)));
187  jk = init_jk[prec];
188  SDL_assert((jk >= 2) && (jk <= 6));
189  jp = jk;
190 
191  /* determine jx,jv,q0, note that 3>q0 */
192  SDL_assert(nx > 0);
193  jx = nx - 1;
194  jv = (e0 - 3) / 24;
195  if (jv < 0)
196  jv = 0;
197  q0 = e0 - 24 * (jv + 1);
198 
199  /* set up f[0] to f[jx+jk] where f[jx+jk] = ipio2[jv+jk] */
200  j = jv - jx;
201  m = jx + jk;
202  for (i = 0; i <= m; i++, j++)
203  f[i] = (j < 0) ? zero : (double) ipio2[j];
204 
205  /* compute q[0],q[1],...q[jk] */
206  for (i = 0; i <= jk; i++) {
207  for (j = 0, fw = 0.0; j <= jx; j++) {
208  const int32_t idx = jx + i - j;
209  SDL_assert(idx >= 0);
210  SDL_assert(idx < 20);
211  SDL_assert(idx <= m);
212  fw += x[j] * f[idx];
213  }
214  q[i] = fw;
215  }
216 
217  jz = jk;
218  recompute:
219  /* distill q[] into iq[] reversingly */
220  for (i = 0, j = jz, z = q[jz]; j > 0; i++, j--) {
221  fw = (double) ((int32_t) (twon24 * z));
222  iq[i] = (int32_t) (z - two24 * fw);
223  z = q[j - 1] + fw;
224  }
225 
226  /* compute n */
227  z = scalbn(z, q0); /* actual value of z */
228  z -= 8.0 * floor(z * 0.125); /* trim off integer >= 8 */
229  n = (int32_t) z;
230  z -= (double) n;
231  ih = 0;
232  if (q0 > 0) { /* need iq[jz-1] to determine n */
233  i = (iq[jz - 1] >> (24 - q0));
234  n += i;
235  iq[jz - 1] -= i << (24 - q0);
236  ih = iq[jz - 1] >> (23 - q0);
237  } else if (q0 == 0)
238  ih = iq[jz - 1] >> 23;
239  else if (z >= 0.5)
240  ih = 2;
241 
242  if (ih > 0) { /* q > 0.5 */
243  n += 1;
244  carry = 0;
245  for (i = 0; i < jz; i++) { /* compute 1-q */
246  j = iq[i];
247  if (carry == 0) {
248  if (j != 0) {
249  carry = 1;
250  iq[i] = 0x1000000 - j;
251  }
252  } else
253  iq[i] = 0xffffff - j;
254  }
255  if (q0 > 0) { /* rare case: chance is 1 in 12 */
256  switch (q0) {
257  case 1:
258  iq[jz - 1] &= 0x7fffff;
259  break;
260  case 2:
261  iq[jz - 1] &= 0x3fffff;
262  break;
263  }
264  }
265  if (ih == 2) {
266  z = one - z;
267  if (carry != 0)
268  z -= scalbn(one, q0);
269  }
270  }
271 
272  /* check if recomputation is needed */
273  if (z == zero) {
274  j = 0;
275  for (i = jz - 1; i >= jk; i--)
276  j |= iq[i];
277  if (j == 0) { /* need recomputation */
278  for (k = 1; iq[jk - k] == 0; k++); /* k = no. of terms needed */
279 
280  for (i = jz + 1; i <= jz + k; i++) { /* add q[jz+1] to q[jz+k] */
281  f[jx + i] = (double) ipio2[jv + i];
282  for (j = 0, fw = 0.0; j <= jx; j++)
283  fw += x[j] * f[jx + i - j];
284  q[i] = fw;
285  }
286  jz += k;
287  goto recompute;
288  }
289  }
290 
291  /* chop off zero terms */
292  if (z == 0.0) {
293  jz -= 1;
294  q0 -= 24;
295  while (iq[jz] == 0) {
296  jz--;
297  q0 -= 24;
298  }
299  } else { /* break z into 24-bit if necessary */
300  z = scalbn(z, -q0);
301  if (z >= two24) {
302  fw = (double) ((int32_t) (twon24 * z));
303  iq[jz] = (int32_t) (z - two24 * fw);
304  jz += 1;
305  q0 += 24;
306  iq[jz] = (int32_t) fw;
307  } else
308  iq[jz] = (int32_t) z;
309  }
310 
311  /* convert integer "bit" chunk to floating-point value */
312  fw = scalbn(one, q0);
313  for (i = jz; i >= 0; i--) {
314  q[i] = fw * (double) iq[i];
315  fw *= twon24;
316  }
317 
318  /* compute PIo2[0,...,jp]*q[jz,...,0] */
319  for (i = jz; i >= 0; i--) {
320  for (fw = 0.0, k = 0; k <= jp && k <= jz - i; k++)
321  fw += PIo2[k] * q[i + k];
322  fq[jz - i] = fw;
323  }
324 
325  /* compress fq[] into y[] */
326  switch (prec) {
327  case 0:
328  fw = 0.0;
329  for (i = jz; i >= 0; i--)
330  fw += fq[i];
331  y[0] = (ih == 0) ? fw : -fw;
332  break;
333  case 1:
334  case 2:
335  fw = 0.0;
336  for (i = jz; i >= 0; i--)
337  fw += fq[i];
338  y[0] = (ih == 0) ? fw : -fw;
339  fw = fq[0] - fw;
340  for (i = 1; i <= jz; i++)
341  fw += fq[i];
342  y[1] = (ih == 0) ? fw : -fw;
343  break;
344  case 3: /* painful */
345  for (i = jz; i > 0; i--) {
346  fw = fq[i - 1] + fq[i];
347  fq[i] += fq[i - 1] - fw;
348  fq[i - 1] = fw;
349  }
350  for (i = jz; i > 1; i--) {
351  fw = fq[i - 1] + fq[i];
352  fq[i] += fq[i - 1] - fw;
353  fq[i - 1] = fw;
354  }
355  for (fw = 0.0, i = jz; i >= 2; i--)
356  fw += fq[i];
357  if (ih == 0) {
358  y[0] = fq[0];
359  y[1] = fq[1];
360  y[2] = fw;
361  } else {
362  y[0] = -fq[0];
363  y[1] = -fq[1];
364  y[2] = -fw;
365  }
366  }
367  return n & 7;
368 }
GLdouble GLdouble z
static double PIo2[]
Definition: k_rem_pio2.c:150
GLint GLint GLint GLint GLint x
Definition: SDL_opengl.h:1574
signed int int32_t
GLdouble GLdouble GLdouble GLdouble q
Definition: SDL_opengl.h:2087
const GLfloat * m
GLfloat f
static double one
Definition: k_rem_pio2.c:167
#define scalbn
Definition: math_private.h:40
static double two24
Definition: k_rem_pio2.c:167
return Display return Display Bool Bool int int int return Display XEvent Bool(*) XPointer return Display return Display Drawable _Xconst char unsigned int unsigned int return Display Pixmap Pixmap XColor XColor unsigned int unsigned int return Display _Xconst char char int char return Display Visual unsigned int int int char unsigned int unsigned int int in j)
Definition: SDL_x11sym.h:50
GLbyte nx
static double twon24
Definition: k_rem_pio2.c:168
GLint GLint GLint GLint GLint GLint y
Definition: SDL_opengl.h:1574
return Display return Display Bool Bool int int int return Display XEvent Bool(*) XPointer return Display return Display Drawable _Xconst char unsigned int unsigned int return Display Pixmap Pixmap XColor XColor unsigned int unsigned int return Display _Xconst char char int char return Display Visual unsigned int int int char unsigned int unsigned int in i)
Definition: SDL_x11sym.h:50
#define SDL_assert(condition)
Definition: SDL_assert.h:169
static double zero
Definition: k_rem_pio2.c:167
GLdouble n
return Display return Display Bool Bool int int int return Display XEvent Bool(*) XPointer return Display return Display Drawable _Xconst char unsigned int unsigned int return Display Pixmap Pixmap XColor XColor unsigned int unsigned int return Display _Xconst char char int char return Display Visual unsigned int int int char unsigned int unsigned int int int return Display Window Cursor return Display Window return Display Drawable GC int int unsigned int unsigned int return Display Drawable GC int int _Xconst char int return Display Drawable GC int int unsigned int unsigned int return Display return Display Cursor return Display GC return XModifierKeymap return char Display Window int return Display return Display Atom return Display Window XWindowAttributes return Display Window return Display XEvent Bool(*) XPointer return Display Window Bool unsigned int int int Window Cursor Time return Display Window int return KeySym return Display _Xconst char Bool return Display _Xconst char return XKeyEvent char int KeySym XComposeStatus return Display int int int XVisualInfo return Display Window int int return _Xconst char return Display XEvent return Display Drawable GC XImage int int int int unsigned int unsigned int return Display Window Window Window int int int int unsigned int return Display Window Window int int return Display Window unsigned int unsigned int return Display Window Bool long XEvent return Display GC unsigned long return Display Window int Time return Display Window Window return Display Window unsigned long return Display Window XSizeHints Display Colormap XColor int return char int XTextProperty return XFontStruct _Xconst char int int int int XCharStruct return Display Window return Display Time return Display Colormap return Display Window Window int int unsigned int unsigned int int int return Display Window int return XExtensionInfo Display char XExtensionHooks int XPointer return XExtensionInfo XExtensionInfo Display return Display return Display unsigned long Display GC Display char long Display xReply int Bool return Display Bool return Display int SDL_X11_XESetEventToWireRetType return Display Window Window Window Window unsigned int return Display XShmSegmentInfo return Display Drawable GC XImage int int int int unsigned int unsigned int Boo k)
Definition: SDL_x11sym.h:213
#define SDL_arraysize(array)
Definition: SDL_stdinc.h:93
#define floor
Definition: math_private.h:37

◆ libm_hidden_proto()

libm_hidden_proto ( scalbn  )

Definition at line 139 of file k_rem_pio2.c.

References PIo2.

142  { 2, 3, 4, 6 }; /* initial value for jk */
143 #else
144  static int init_jk[] = { 2, 3, 4, 6 };

Variable Documentation

◆ one

double one = 1.0
static

Definition at line 167 of file k_rem_pio2.c.

Referenced by __kernel_rem_pio2().

◆ PIo2

double PIo2[]
static
Initial value:
= {
1.57079625129699707031e+00,
7.54978941586159635335e-08,
5.39030252995776476554e-15,
3.28200341580791294123e-22,
1.27065575308067607349e-29,
1.22933308981111328932e-36,
2.73370053816464559624e-44,
2.16741683877804819444e-51,
}

Definition at line 150 of file k_rem_pio2.c.

Referenced by libm_hidden_proto().

◆ two24

double two24 = 1.67772160000000000000e+07
static

Definition at line 167 of file k_rem_pio2.c.

Referenced by __kernel_rem_pio2().

◆ twon24

double twon24 = 5.96046447753906250000e-08
static

Definition at line 168 of file k_rem_pio2.c.

Referenced by __kernel_rem_pio2().

◆ zero

double zero = 0.0
static

Definition at line 167 of file k_rem_pio2.c.

Referenced by __kernel_rem_pio2().