/[jscoverage]/trunk/js/src/fdlibm/e_fmod.c
ViewVC logotype

Contents of /trunk/js/src/fdlibm/e_fmod.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2 - (show annotations)
Wed Aug 1 13:51:53 2007 UTC (14 years, 9 months ago) by siliconforks
File MIME type: text/plain
File size: 5458 byte(s)
Initial import.

1 /* -*- Mode: C; tab-width: 8; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 *
3 * ***** BEGIN LICENSE BLOCK *****
4 * Version: MPL 1.1/GPL 2.0/LGPL 2.1
5 *
6 * The contents of this file are subject to the Mozilla Public License Version
7 * 1.1 (the "License"); you may not use this file except in compliance with
8 * the License. You may obtain a copy of the License at
9 * http://www.mozilla.org/MPL/
10 *
11 * Software distributed under the License is distributed on an "AS IS" basis,
12 * WITHOUT WARRANTY OF ANY KIND, either express or implied. See the License
13 * for the specific language governing rights and limitations under the
14 * License.
15 *
16 * The Original Code is Mozilla Communicator client code, released
17 * March 31, 1998.
18 *
19 * The Initial Developer of the Original Code is
20 * Sun Microsystems, Inc.
21 * Portions created by the Initial Developer are Copyright (C) 1998
22 * the Initial Developer. All Rights Reserved.
23 *
24 * Contributor(s):
25 *
26 * Alternatively, the contents of this file may be used under the terms of
27 * either of the GNU General Public License Version 2 or later (the "GPL"),
28 * or the GNU Lesser General Public License Version 2.1 or later (the "LGPL"),
29 * in which case the provisions of the GPL or the LGPL are applicable instead
30 * of those above. If you wish to allow use of your version of this file only
31 * under the terms of either the GPL or the LGPL, and not to allow others to
32 * use your version of this file under the terms of the MPL, indicate your
33 * decision by deleting the provisions above and replace them with the notice
34 * and other provisions required by the GPL or the LGPL. If you do not delete
35 * the provisions above, a recipient may use your version of this file under
36 * the terms of any one of the MPL, the GPL or the LGPL.
37 *
38 * ***** END LICENSE BLOCK ***** */
39
40 /* @(#)e_fmod.c 1.3 95/01/18 */
41 /*
42 * ====================================================
43 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
44 *
45 * Developed at SunSoft, a Sun Microsystems, Inc. business.
46 * Permission to use, copy, modify, and distribute this
47 * software is freely granted, provided that this notice
48 * is preserved.
49 * ====================================================
50 */
51
52 /*
53 * __ieee754_fmod(x,y)
54 * Return x mod y in exact arithmetic
55 * Method: shift and subtract
56 */
57
58 #include "fdlibm.h"
59
60 #ifdef __STDC__
61 static const double one = 1.0, Zero[] = {0.0, -0.0,};
62 #else
63 static double one = 1.0, Zero[] = {0.0, -0.0,};
64 #endif
65
66 #ifdef __STDC__
67 double __ieee754_fmod(double x, double y)
68 #else
69 double __ieee754_fmod(x,y)
70 double x,y ;
71 #endif
72 {
73 fd_twoints ux, uy;
74 int n,hx,hy,hz,ix,iy,sx,i;
75 unsigned lx,ly,lz;
76
77 ux.d = x; uy.d = y;
78 hx = __HI(ux); /* high word of x */
79 lx = __LO(ux); /* low word of x */
80 hy = __HI(uy); /* high word of y */
81 ly = __LO(uy); /* low word of y */
82 sx = hx&0x80000000; /* sign of x */
83 hx ^=sx; /* |x| */
84 hy &= 0x7fffffff; /* |y| */
85
86 /* purge off exception values */
87 if((hy|ly)==0||(hx>=0x7ff00000)|| /* y=0,or x not finite */
88 ((hy|((ly|-(int)ly)>>31))>0x7ff00000)) /* or y is NaN */
89 return (x*y)/(x*y);
90 if(hx<=hy) {
91 if((hx<hy)||(lx<ly)) return x; /* |x|<|y| return x */
92 if(lx==ly)
93 return Zero[(unsigned)sx>>31]; /* |x|=|y| return x*0*/
94 }
95
96 /* determine ix = ilogb(x) */
97 if(hx<0x00100000) { /* subnormal x */
98 if(hx==0) {
99 for (ix = -1043, i=lx; i>0; i<<=1) ix -=1;
100 } else {
101 for (ix = -1022,i=(hx<<11); i>0; i<<=1) ix -=1;
102 }
103 } else ix = (hx>>20)-1023;
104
105 /* determine iy = ilogb(y) */
106 if(hy<0x00100000) { /* subnormal y */
107 if(hy==0) {
108 for (iy = -1043, i=ly; i>0; i<<=1) iy -=1;
109 } else {
110 for (iy = -1022,i=(hy<<11); i>0; i<<=1) iy -=1;
111 }
112 } else iy = (hy>>20)-1023;
113
114 /* set up {hx,lx}, {hy,ly} and align y to x */
115 if(ix >= -1022)
116 hx = 0x00100000|(0x000fffff&hx);
117 else { /* subnormal x, shift x to normal */
118 n = -1022-ix;
119 if(n<=31) {
120 hx = (hx<<n)|(lx>>(32-n));
121 lx <<= n;
122 } else {
123 hx = lx<<(n-32);
124 lx = 0;
125 }
126 }
127 if(iy >= -1022)
128 hy = 0x00100000|(0x000fffff&hy);
129 else { /* subnormal y, shift y to normal */
130 n = -1022-iy;
131 if(n<=31) {
132 hy = (hy<<n)|(ly>>(32-n));
133 ly <<= n;
134 } else {
135 hy = ly<<(n-32);
136 ly = 0;
137 }
138 }
139
140 /* fix point fmod */
141 n = ix - iy;
142 while(n--) {
143 hz=hx-hy;lz=lx-ly; if(lx<ly) hz -= 1;
144 if(hz<0){hx = hx+hx+(lx>>31); lx = lx+lx;}
145 else {
146 if((hz|lz)==0) /* return sign(x)*0 */
147 return Zero[(unsigned)sx>>31];
148 hx = hz+hz+(lz>>31); lx = lz+lz;
149 }
150 }
151 hz=hx-hy;lz=lx-ly; if(lx<ly) hz -= 1;
152 if(hz>=0) {hx=hz;lx=lz;}
153
154 /* convert back to floating value and restore the sign */
155 if((hx|lx)==0) /* return sign(x)*0 */
156 return Zero[(unsigned)sx>>31];
157 while(hx<0x00100000) { /* normalize x */
158 hx = hx+hx+(lx>>31); lx = lx+lx;
159 iy -= 1;
160 }
161 if(iy>= -1022) { /* normalize output */
162 hx = ((hx-0x00100000)|((iy+1023)<<20));
163 ux.d = x;
164 __HI(ux) = hx|sx;
165 __LO(ux) = lx;
166 x = ux.d;
167 } else { /* subnormal output */
168 n = -1022 - iy;
169 if(n<=20) {
170 lx = (lx>>n)|((unsigned)hx<<(32-n));
171 hx >>= n;
172 } else if (n<=31) {
173 lx = (hx<<(32-n))|(lx>>n); hx = sx;
174 } else {
175 lx = hx>>(n-32); hx = sx;
176 }
177 ux.d = x;
178 __HI(ux) = hx|sx;
179 __LO(ux) = lx;
180 x = ux.d;
181 x *= one; /* create necessary signal */
182 }
183 return x; /* exact output */
184 }

  ViewVC Help
Powered by ViewVC 1.1.24