|
|
1.1 root 1: /*
2: * Copyright (c) 1989 The Regents of the University of California.
3: * All rights reserved.
4: *
5: * Redistribution and use in source and binary forms are permitted
6: * provided that: (1) source distributions retain this entire copyright
7: * notice and comment, and (2) distributions including binaries display
8: * the following acknowledgement: ``This product includes software
9: * developed by the University of California, Berkeley and its contributors''
10: * in the documentation or other materials provided with the distribution
11: * and in all advertising materials mentioning features or use of this
12: * software. Neither the name of the University nor the names of its
13: * contributors may be used to endorse or promote products derived
14: * from this software without specific prior written permission.
15: * THIS SOFTWARE IS PROVIDED ``AS IS'' AND WITHOUT ANY EXPRESS OR
16: * IMPLIED WARRANTIES, INCLUDING, WITHOUT LIMITATION, THE IMPLIED
17: * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
18: */
19:
20: #ifndef lint
21: static char sccsid[] = "@(#)fmod.c 5.2 (Berkeley) 6/1/90";
22: #endif /* not lint */
23:
24: /* fmod.c
25: *
26: * SYNOPSIS
27: *
28: * #include <math.h>
29: * double fmod(double x, double y)
30: *
31: * DESCRIPTION
32: *
33: * The fmod function computes the floating-point remainder of x/y.
34: *
35: * RETURNS
36: *
37: * The fmod function returns the value x-i*y, for some integer i
38: * such that, if y is nonzero, the result has the same sign as x and
39: * magnitude less than the magnitude of y.
40: *
41: * On a VAX or CCI,
42: *
43: * fmod(x,0) traps/faults on floating-point divided-by-zero.
44: *
45: * On IEEE-754 conforming machines with "isnan()" primitive,
46: *
47: * fmod(x,0), fmod(INF,y) are invalid operations and NaN is returned.
48: *
49: */
50: #if !defined(vax) && !defined(tahoe)
51: extern int isnan(),finite();
52: #endif /* !defined(vax) && !defined(tahoe) */
53: extern double frexp(),ldexp(),fabs();
54:
55: #ifdef TEST_FMOD
56: static double
57: _fmod(x,y)
58: #else /* TEST_FMOD */
59: double
60: fmod(x,y)
61: #endif /* TEST_FMOD */
62: double x,y;
63: {
64: int ir,iy;
65: double r,w;
66:
67: if (y == (double)0
68: #if !defined(vax) && !defined(tahoe) /* per "fmod" manual entry, SunOS 4.0 */
69: || isnan(y) || !finite(x)
70: #endif /* !defined(vax) && !defined(tahoe) */
71: )
72: return (x*y)/(x*y);
73:
74: r = fabs(x);
75: y = fabs(y);
76: (void)frexp(y,&iy);
77: while (r >= y) {
78: (void)frexp(r,&ir);
79: w = ldexp(y,ir-iy);
80: r -= w <= r ? w : w*(double)0.5;
81: }
82: return x >= (double)0 ? r : -r;
83: }
84:
85: #ifdef TEST_FMOD
86: extern long random();
87: extern double fmod();
88:
89: #define NTEST 10000
90: #define NCASES 3
91:
92: static int nfail = 0;
93:
94: static void
95: doit(x,y)
96: double x,y;
97: {
98: double ro = fmod(x,y),rn = _fmod(x,y);
99: if (ro != rn) {
100: (void)printf(" x = 0x%08.8x %08.8x (%24.16e)\n",x,x);
101: (void)printf(" y = 0x%08.8x %08.8x (%24.16e)\n",y,y);
102: (void)printf(" fmod = 0x%08.8x %08.8x (%24.16e)\n",ro,ro);
103: (void)printf("_fmod = 0x%08.8x %08.8x (%24.16e)\n",rn,rn);
104: (void)printf("\n");
105: }
106: }
107:
108: main()
109: {
110: register int i,cases;
111: double x,y;
112:
113: srandom(12345);
114: for (i = 0; i < NTEST; i++) {
115: x = (double)random();
116: y = (double)random();
117: for (cases = 0; cases < NCASES; cases++) {
118: switch (cases) {
119: case 0:
120: break;
121: case 1:
122: y = (double)1/y; break;
123: case 2:
124: x = (double)1/x; break;
125: default:
126: abort(); break;
127: }
128: doit(x,y);
129: doit(x,-y);
130: doit(-x,y);
131: doit(-x,-y);
132: }
133: }
134: if (nfail)
135: (void)printf("Number of failures: %d (out of a total of %d)\n",
136: nfail,NTEST*NCASES*4);
137: else
138: (void)printf("No discrepancies were found\n");
139: exit(0);
140: }
141: #endif /* TEST_FMOD */
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.