1 /* IEEE754 floating point arithmetic 2 * single precision 3 */ 4 /* 5 * MIPS floating point support 6 * Copyright (C) 1994-2000 Algorithmics Ltd. 7 * 8 * ######################################################################## 9 * 10 * This program is free software; you can distribute it and/or modify it 11 * under the terms of the GNU General Public License (Version 2) as 12 * published by the Free Software Foundation. 13 * 14 * This program is distributed in the hope it will be useful, but WITHOUT 15 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 16 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 17 * for more details. 18 * 19 * You should have received a copy of the GNU General Public License along 20 * with this program; if not, write to the Free Software Foundation, Inc., 21 * 59 Temple Place - Suite 330, Boston MA 02111-1307, USA. 22 * 23 * ######################################################################## 24 */ 25 26 27 #include "ieee754sp.h" 28 29 int ieee754sp_class(ieee754sp x) 30 { 31 COMPXSP; 32 EXPLODEXSP; 33 return xc; 34 } 35 36 int ieee754sp_isnan(ieee754sp x) 37 { 38 return ieee754sp_class(x) >= IEEE754_CLASS_SNAN; 39 } 40 41 int ieee754sp_issnan(ieee754sp x) 42 { 43 assert(ieee754sp_isnan(x)); 44 return (SPMANT(x) & SP_MBIT(SP_MBITS-1)); 45 } 46 47 48 ieee754sp ieee754sp_xcpt(ieee754sp r, const char *op, ...) 49 { 50 struct ieee754xctx ax; 51 52 if (!TSTX()) 53 return r; 54 55 ax.op = op; 56 ax.rt = IEEE754_RT_SP; 57 ax.rv.sp = r; 58 va_start(ax.ap, op); 59 ieee754_xcpt(&ax); 60 va_end(ax.ap); 61 return ax.rv.sp; 62 } 63 64 ieee754sp ieee754sp_nanxcpt(ieee754sp r, const char *op, ...) 65 { 66 struct ieee754xctx ax; 67 68 assert(ieee754sp_isnan(r)); 69 70 if (!ieee754sp_issnan(r)) /* QNAN does not cause invalid op !! */ 71 return r; 72 73 if (!SETANDTESTCX(IEEE754_INVALID_OPERATION)) { 74 /* not enabled convert to a quiet NaN */ 75 SPMANT(r) &= (~SP_MBIT(SP_MBITS-1)); 76 if (ieee754sp_isnan(r)) 77 return r; 78 else 79 return ieee754sp_indef(); 80 } 81 82 ax.op = op; 83 ax.rt = 0; 84 ax.rv.sp = r; 85 va_start(ax.ap, op); 86 ieee754_xcpt(&ax); 87 va_end(ax.ap); 88 return ax.rv.sp; 89 } 90 91 ieee754sp ieee754sp_bestnan(ieee754sp x, ieee754sp y) 92 { 93 assert(ieee754sp_isnan(x)); 94 assert(ieee754sp_isnan(y)); 95 96 if (SPMANT(x) > SPMANT(y)) 97 return x; 98 else 99 return y; 100 } 101 102 103 static unsigned get_rounding(int sn, unsigned xm) 104 { 105 /* inexact must round of 3 bits 106 */ 107 if (xm & (SP_MBIT(3) - 1)) { 108 switch (ieee754_csr.rm) { 109 case IEEE754_RZ: 110 break; 111 case IEEE754_RN: 112 xm += 0x3 + ((xm >> 3) & 1); 113 /* xm += (xm&0x8)?0x4:0x3 */ 114 break; 115 case IEEE754_RU: /* toward +Infinity */ 116 if (!sn) /* ?? */ 117 xm += 0x8; 118 break; 119 case IEEE754_RD: /* toward -Infinity */ 120 if (sn) /* ?? */ 121 xm += 0x8; 122 break; 123 } 124 } 125 return xm; 126 } 127 128 129 /* generate a normal/denormal number with over,under handling 130 * sn is sign 131 * xe is an unbiased exponent 132 * xm is 3bit extended precision value. 133 */ 134 ieee754sp ieee754sp_format(int sn, int xe, unsigned xm) 135 { 136 assert(xm); /* we don't gen exact zeros (probably should) */ 137 138 assert((xm >> (SP_MBITS + 1 + 3)) == 0); /* no execess */ 139 assert(xm & (SP_HIDDEN_BIT << 3)); 140 141 if (xe < SP_EMIN) { 142 /* strip lower bits */ 143 int es = SP_EMIN - xe; 144 145 if (ieee754_csr.nod) { 146 SETCX(IEEE754_UNDERFLOW); 147 SETCX(IEEE754_INEXACT); 148 149 switch(ieee754_csr.rm) { 150 case IEEE754_RN: 151 case IEEE754_RZ: 152 return ieee754sp_zero(sn); 153 case IEEE754_RU: /* toward +Infinity */ 154 if(sn == 0) 155 return ieee754sp_min(0); 156 else 157 return ieee754sp_zero(1); 158 case IEEE754_RD: /* toward -Infinity */ 159 if(sn == 0) 160 return ieee754sp_zero(0); 161 else 162 return ieee754sp_min(1); 163 } 164 } 165 166 if (xe == SP_EMIN - 1 167 && get_rounding(sn, xm) >> (SP_MBITS + 1 + 3)) 168 { 169 /* Not tiny after rounding */ 170 SETCX(IEEE754_INEXACT); 171 xm = get_rounding(sn, xm); 172 xm >>= 1; 173 /* Clear grs bits */ 174 xm &= ~(SP_MBIT(3) - 1); 175 xe++; 176 } 177 else { 178 /* sticky right shift es bits 179 */ 180 SPXSRSXn(es); 181 assert((xm & (SP_HIDDEN_BIT << 3)) == 0); 182 assert(xe == SP_EMIN); 183 } 184 } 185 if (xm & (SP_MBIT(3) - 1)) { 186 SETCX(IEEE754_INEXACT); 187 if ((xm & (SP_HIDDEN_BIT << 3)) == 0) { 188 SETCX(IEEE754_UNDERFLOW); 189 } 190 191 /* inexact must round of 3 bits 192 */ 193 xm = get_rounding(sn, xm); 194 /* adjust exponent for rounding add overflowing 195 */ 196 if (xm >> (SP_MBITS + 1 + 3)) { 197 /* add causes mantissa overflow */ 198 xm >>= 1; 199 xe++; 200 } 201 } 202 /* strip grs bits */ 203 xm >>= 3; 204 205 assert((xm >> (SP_MBITS + 1)) == 0); /* no execess */ 206 assert(xe >= SP_EMIN); 207 208 if (xe > SP_EMAX) { 209 SETCX(IEEE754_OVERFLOW); 210 SETCX(IEEE754_INEXACT); 211 /* -O can be table indexed by (rm,sn) */ 212 switch (ieee754_csr.rm) { 213 case IEEE754_RN: 214 return ieee754sp_inf(sn); 215 case IEEE754_RZ: 216 return ieee754sp_max(sn); 217 case IEEE754_RU: /* toward +Infinity */ 218 if (sn == 0) 219 return ieee754sp_inf(0); 220 else 221 return ieee754sp_max(1); 222 case IEEE754_RD: /* toward -Infinity */ 223 if (sn == 0) 224 return ieee754sp_max(0); 225 else 226 return ieee754sp_inf(1); 227 } 228 } 229 /* gen norm/denorm/zero */ 230 231 if ((xm & SP_HIDDEN_BIT) == 0) { 232 /* we underflow (tiny/zero) */ 233 assert(xe == SP_EMIN); 234 if (ieee754_csr.mx & IEEE754_UNDERFLOW) 235 SETCX(IEEE754_UNDERFLOW); 236 return buildsp(sn, SP_EMIN - 1 + SP_EBIAS, xm); 237 } else { 238 assert((xm >> (SP_MBITS + 1)) == 0); /* no execess */ 239 assert(xm & SP_HIDDEN_BIT); 240 241 return buildsp(sn, xe + SP_EBIAS, xm & ~SP_HIDDEN_BIT); 242 } 243 } 244