1 /* 2 * CDDL HEADER START 3 * 4 * The contents of this file are subject to the terms of the 5 * Common Development and Distribution License, Version 1.0 only 6 * (the "License"). You may not use this file except in compliance 7 * with the License. 8 * 9 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE 10 * or http://www.opensolaris.org/os/licensing. 11 * See the License for the specific language governing permissions 12 * and limitations under the License. 13 * 14 * When distributing Covered Code, include this CDDL HEADER in each 15 * file and include the License file at usr/src/OPENSOLARIS.LICENSE. 16 * If applicable, add the following below this CDDL HEADER, with the 17 * fields enclosed by brackets "[]" replaced with your own identifying 18 * information: Portions Copyright [yyyy] [name of copyright owner] 19 * 20 * CDDL HEADER END 21 */ 22 /* 23 * Copyright (c) 1988 by Sun Microsystems, Inc. 24 */ 25 26 #ident "%Z%%M% %I% %E% SMI" /* SunOS-4.1 1.9 88/11/30 */ 27 28 #include <sys/fpu/fpu_simulator.h> 29 #include <sys/fpu/globals.h> 30 31 void 32 _fp_div(pfpsd, px, py, pz) 33 fp_simd_type *pfpsd; 34 unpacked *px, *py, *pz; 35 { 36 unsigned r[4], *y, q, c; 37 int n; 38 39 *pz = *px; 40 41 if ((py->fpclass >= fp_quiet) || (px->fpclass >= fp_quiet)) { 42 if (py->fpclass >= px->fpclass) *pz = *py; 43 return; 44 } 45 46 pz->sign = px->sign ^ py->sign; 47 switch (px->fpclass) { 48 case fp_quiet: 49 case fp_signaling: 50 return; 51 case fp_zero: 52 case fp_infinity: 53 if (px->fpclass == py->fpclass) { /* 0/0 or inf/inf */ 54 fpu_error_nan(pfpsd, pz); 55 pz->fpclass = fp_quiet; 56 } 57 return; 58 case fp_normal: 59 switch (py->fpclass) { 60 case fp_zero: /* number/0 */ 61 fpu_set_exception(pfpsd, fp_division); 62 pz->fpclass = fp_infinity; 63 return; 64 case fp_infinity: /* number/inf */ 65 pz->fpclass = fp_zero; 66 return; 67 } 68 } 69 70 /* Now x and y are both normal or subnormal. */ 71 72 r[0] = px->significand[0]; 73 r[1] = px->significand[1]; 74 r[2] = px->significand[2]; 75 r[3] = px->significand[3]; 76 y = py->significand; 77 78 if (fpu_cmpli(r, y, 4) >= 0) 79 pz->exponent = px->exponent - py->exponent; 80 else 81 pz->exponent = px->exponent - py->exponent - 1; 82 83 q = 0; 84 while (q < 0x10000) { /* generate quo[0] */ 85 q <<= 1; 86 if (fpu_cmpli(r, y, 4) >= 0) { 87 q += 1; /* if r>y do r-=y and q+=1 */ 88 c = 0; 89 c = fpu_sub3wc(&r[3], r[3], y[3], c); 90 c = fpu_sub3wc(&r[2], r[2], y[2], c); 91 c = fpu_sub3wc(&r[1], r[1], y[1], c); 92 c = fpu_sub3wc(&r[0], r[0], y[0], c); 93 } 94 r[0] = (r[0]<<1)|((r[1]&0x80000000)>>31); /* r << 1 */ 95 r[1] = (r[1]<<1)|((r[2]&0x80000000)>>31); 96 r[2] = (r[2]<<1)|((r[3]&0x80000000)>>31); 97 r[3] = (r[3]<<1); 98 } 99 pz->significand[0] = q; 100 q = 0; /* generate quo[1] */ 101 n = 32; 102 while (n--) { 103 q <<= 1; 104 if (fpu_cmpli(r, y, 4) >= 0) { 105 q += 1; /* if r>y do r-=y and q+=1 */ 106 c = 0; 107 c = fpu_sub3wc(&r[3], r[3], y[3], c); 108 c = fpu_sub3wc(&r[2], r[2], y[2], c); 109 c = fpu_sub3wc(&r[1], r[1], y[1], c); 110 c = fpu_sub3wc(&r[0], r[0], y[0], c); 111 } 112 r[0] = (r[0]<<1)|((r[1]&0x80000000)>>31); /* r << 1 */ 113 r[1] = (r[1]<<1)|((r[2]&0x80000000)>>31); 114 r[2] = (r[2]<<1)|((r[3]&0x80000000)>>31); 115 r[3] = (r[3]<<1); 116 } 117 pz->significand[1] = q; 118 q = 0; /* generate quo[2] */ 119 n = 32; 120 while (n--) { 121 q <<= 1; 122 if (fpu_cmpli(r, y, 4) >= 0) { 123 q += 1; /* if r>y do r-=y and q+=1 */ 124 c = 0; 125 c = fpu_sub3wc(&r[3], r[3], y[3], c); 126 c = fpu_sub3wc(&r[2], r[2], y[2], c); 127 c = fpu_sub3wc(&r[1], r[1], y[1], c); 128 c = fpu_sub3wc(&r[0], r[0], y[0], c); 129 } 130 r[0] = (r[0]<<1)|((r[1]&0x80000000)>>31); /* r << 1 */ 131 r[1] = (r[1]<<1)|((r[2]&0x80000000)>>31); 132 r[2] = (r[2]<<1)|((r[3]&0x80000000)>>31); 133 r[3] = (r[3]<<1); 134 } 135 pz->significand[2] = q; 136 q = 0; /* generate quo[3] */ 137 n = 32; 138 while (n--) { 139 q <<= 1; 140 if (fpu_cmpli(r, y, 4) >= 0) { 141 q += 1; /* if r>y do r-=y and q+=1 */ 142 c = 0; 143 c = fpu_sub3wc(&r[3], r[3], y[3], c); 144 c = fpu_sub3wc(&r[2], r[2], y[2], c); 145 c = fpu_sub3wc(&r[1], r[1], y[1], c); 146 c = fpu_sub3wc(&r[0], r[0], y[0], c); 147 } 148 r[0] = (r[0]<<1)|((r[1]&0x80000000)>>31); /* r << 1 */ 149 r[1] = (r[1]<<1)|((r[2]&0x80000000)>>31); 150 r[2] = (r[2]<<1)|((r[3]&0x80000000)>>31); 151 r[3] = (r[3]<<1); 152 } 153 pz->significand[3] = q; 154 if ((r[0]|r[1]|r[2]|r[3]) == 0) pz->sticky = pz->rounded = 0; 155 else { 156 pz->sticky = 1; /* half way case won't occur */ 157 if (fpu_cmpli(r, y, 4) >= 0) pz->rounded = 1; 158 } 159 } 160 161 void 162 _fp_sqrt(pfpsd, px, pz) 163 fp_simd_type *pfpsd; 164 unpacked *px, *pz; 165 { /* *pz gets sqrt(*px) */ 166 167 unsigned *x, r, c, q, t[4], s[4]; 168 *pz = *px; 169 switch (px->fpclass) { 170 case fp_quiet: 171 case fp_signaling: 172 case fp_zero: 173 return; 174 case fp_infinity: 175 if (px->sign == 1) { /* sqrt(-inf) */ 176 fpu_error_nan(pfpsd, pz); 177 pz->fpclass = fp_quiet; 178 } 179 return; 180 case fp_normal: 181 if (px->sign == 1) { /* sqrt(-norm) */ 182 fpu_error_nan(pfpsd, pz); 183 pz->fpclass = fp_quiet; 184 return; 185 } 186 } 187 188 /* Now x is normal. */ 189 x = px->significand; 190 if (px->exponent & 1) { 191 /* 192 * sqrt(1.f * 2**odd) = sqrt (2.+2f) 193 * 2**(odd-1)/2 194 */ 195 pz->exponent = (px->exponent - 1) / 2; 196 x[0] = (x[0]<<1)|((x[1]&0x80000000)>>31); /* x<<1 */ 197 x[1] = (x[1]<<1)|((x[2]&0x80000000)>>31); 198 x[2] = (x[2]<<1)|((x[3]&0x80000000)>>31); 199 x[3] = (x[3]<<1); 200 } else { 201 /* 202 * sqrt(1.f * 2**even) = sqrt (1.f) 203 * 2**(even)/2 204 */ 205 pz->exponent = px->exponent / 2; 206 } 207 s[0] = s[1] = s[2] = s[3] = t[0] = t[1] = t[2] = t[3] = 0; 208 q = 0; 209 r = 0x00010000; 210 while (r != 0) { /* compute sqrt[0] */ 211 t[0] = s[0] + r; 212 if (t[0] <= x[0]) { 213 s[0] = t[0] + r; 214 x[0] -= t[0]; 215 q += r; 216 } 217 x[0] = (x[0]<<1)|((x[1]&0x80000000)>>31); /* x<<1 */ 218 x[1] = (x[1]<<1)|((x[2]&0x80000000)>>31); 219 x[2] = (x[2]<<1)|((x[3]&0x80000000)>>31); 220 x[3] = (x[3]<<1); 221 r >>= 1; 222 } 223 pz->significand[0] = q; 224 q = 0; 225 r = (unsigned)0x80000000; 226 while (r != 0) { /* compute sqrt[1] */ 227 t[1] = s[1] + r; /* no carry */ 228 t[0] = s[0]; 229 if (fpu_cmpli(t, x, 2) <= 0) { 230 c = 0; 231 c = fpu_add3wc(&s[1], t[1], r, c); 232 c = fpu_add3wc(&s[0], t[0], 0, c); 233 c = 0; 234 c = fpu_sub3wc(&x[1], x[1], t[1], c); 235 c = fpu_sub3wc(&x[0], x[0], t[0], c); 236 q += r; 237 } 238 x[0] = (x[0]<<1)|((x[1]&0x80000000)>>31); /* x<<1 */ 239 x[1] = (x[1]<<1)|((x[2]&0x80000000)>>31); 240 x[2] = (x[2]<<1)|((x[3]&0x80000000)>>31); 241 x[3] = (x[3]<<1); 242 r >>= 1; 243 } 244 pz->significand[1] = q; 245 q = 0; 246 r = (unsigned)0x80000000; 247 while (r != 0) { /* compute sqrt[2] */ 248 t[2] = s[2] + r; /* no carry */ 249 t[1] = s[1]; 250 t[0] = s[0]; 251 if (fpu_cmpli(t, x, 3) <= 0) { 252 c = 0; 253 c = fpu_add3wc(&s[2], t[2], r, c); 254 c = fpu_add3wc(&s[1], t[1], 0, c); 255 c = fpu_add3wc(&s[0], t[0], 0, c); 256 c = 0; 257 c = fpu_sub3wc(&x[2], x[2], t[2], c); 258 c = fpu_sub3wc(&x[1], x[1], t[1], c); 259 c = fpu_sub3wc(&x[0], x[0], t[0], c); 260 q += r; 261 } 262 x[0] = (x[0]<<1)|((x[1]&0x80000000)>>31); /* x<<1 */ 263 x[1] = (x[1]<<1)|((x[2]&0x80000000)>>31); 264 x[2] = (x[2]<<1)|((x[3]&0x80000000)>>31); 265 x[3] = (x[3]<<1); 266 r >>= 1; 267 } 268 pz->significand[2] = q; 269 q = 0; 270 r = (unsigned)0x80000000; 271 while (r != 0) { /* compute sqrt[3] */ 272 t[3] = s[3] + r; /* no carry */ 273 t[2] = s[2]; 274 t[1] = s[1]; 275 t[0] = s[0]; 276 if (fpu_cmpli(t, x, 4) <= 0) { 277 c = 0; 278 c = fpu_add3wc(&s[3], t[3], r, c); 279 c = fpu_add3wc(&s[2], t[2], 0, c); 280 c = fpu_add3wc(&s[1], t[1], 0, c); 281 c = fpu_add3wc(&s[0], t[0], 0, c); 282 c = 0; 283 c = fpu_sub3wc(&x[3], x[3], t[3], c); 284 c = fpu_sub3wc(&x[2], x[2], t[2], c); 285 c = fpu_sub3wc(&x[1], x[1], t[1], c); 286 c = fpu_sub3wc(&x[0], x[0], t[0], c); 287 q += r; 288 } 289 x[0] = (x[0]<<1)|((x[1]&0x80000000)>>31); /* x<<1 */ 290 x[1] = (x[1]<<1)|((x[2]&0x80000000)>>31); 291 x[2] = (x[2]<<1)|((x[3]&0x80000000)>>31); 292 x[3] = (x[3]<<1); 293 r >>= 1; 294 } 295 pz->significand[3] = q; 296 if ((x[0]|x[1]|x[2]|x[3]) == 0) { 297 pz->sticky = pz->rounded = 0; 298 } else { 299 pz->sticky = 1; 300 if (fpu_cmpli(s, x, 4) < 0) 301 pz->rounded = 1; 302 else 303 pz->rounded = 0; 304 } 305 } 306