1 /* 2 * Copyright (c) 1989, 1993 3 * The Regents of the University of California. All rights reserved. 4 * 5 * %sccs.include.redist.c% 6 */ 7 8 #ifndef lint 9 static char sccsid[] = "@(#)fmod.c 8.1 (Berkeley) 06/04/93"; 10 #endif /* not lint */ 11 12 /* fmod.c 13 * 14 * SYNOPSIS 15 * 16 * #include <math.h> 17 * double fmod(double x, double y) 18 * 19 * DESCRIPTION 20 * 21 * The fmod function computes the floating-point remainder of x/y. 22 * 23 * RETURNS 24 * 25 * The fmod function returns the value x-i*y, for some integer i 26 * such that, if y is nonzero, the result has the same sign as x and 27 * magnitude less than the magnitude of y. 28 * 29 * On a VAX or CCI, 30 * 31 * fmod(x,0) traps/faults on floating-point divided-by-zero. 32 * 33 * On IEEE-754 conforming machines with "isnan()" primitive, 34 * 35 * fmod(x,0), fmod(INF,y) are invalid operations and NaN is returned. 36 * 37 */ 38 #if !defined(vax) && !defined(tahoe) 39 extern int isnan(),finite(); 40 #endif /* !defined(vax) && !defined(tahoe) */ 41 extern double frexp(),ldexp(),fabs(); 42 43 #ifdef TEST_FMOD 44 static double 45 _fmod(x,y) 46 #else /* TEST_FMOD */ 47 double 48 fmod(x,y) 49 #endif /* TEST_FMOD */ 50 double x,y; 51 { 52 int ir,iy; 53 double r,w; 54 55 if (y == (double)0 56 #if !defined(vax) && !defined(tahoe) /* per "fmod" manual entry, SunOS 4.0 */ 57 || isnan(y) || !finite(x) 58 #endif /* !defined(vax) && !defined(tahoe) */ 59 ) 60 return (x*y)/(x*y); 61 62 r = fabs(x); 63 y = fabs(y); 64 (void)frexp(y,&iy); 65 while (r >= y) { 66 (void)frexp(r,&ir); 67 w = ldexp(y,ir-iy); 68 r -= w <= r ? w : w*(double)0.5; 69 } 70 return x >= (double)0 ? r : -r; 71 } 72 73 #ifdef TEST_FMOD 74 extern long random(); 75 extern double fmod(); 76 77 #define NTEST 10000 78 #define NCASES 3 79 80 static int nfail = 0; 81 82 static void 83 doit(x,y) 84 double x,y; 85 { 86 double ro = fmod(x,y),rn = _fmod(x,y); 87 if (ro != rn) { 88 (void)printf(" x = 0x%08.8x %08.8x (%24.16e)\n",x,x); 89 (void)printf(" y = 0x%08.8x %08.8x (%24.16e)\n",y,y); 90 (void)printf(" fmod = 0x%08.8x %08.8x (%24.16e)\n",ro,ro); 91 (void)printf("_fmod = 0x%08.8x %08.8x (%24.16e)\n",rn,rn); 92 (void)printf("\n"); 93 } 94 } 95 96 main() 97 { 98 register int i,cases; 99 double x,y; 100 101 srandom(12345); 102 for (i = 0; i < NTEST; i++) { 103 x = (double)random(); 104 y = (double)random(); 105 for (cases = 0; cases < NCASES; cases++) { 106 switch (cases) { 107 case 0: 108 break; 109 case 1: 110 y = (double)1/y; break; 111 case 2: 112 x = (double)1/x; break; 113 default: 114 abort(); break; 115 } 116 doit(x,y); 117 doit(x,-y); 118 doit(-x,y); 119 doit(-x,-y); 120 } 121 } 122 if (nfail) 123 (void)printf("Number of failures: %d (out of a total of %d)\n", 124 nfail,NTEST*NCASES*4); 125 else 126 (void)printf("No discrepancies were found\n"); 127 exit(0); 128 } 129 #endif /* TEST_FMOD */ 130