1 /* IEEE754 floating point arithmetic 2 * double precision: common utilities 3 */ 4 /* 5 * MIPS floating point support 6 * Copyright (C) 1994-2000 Algorithmics Ltd. 7 * 8 * This program is free software; you can distribute it and/or modify it 9 * under the terms of the GNU General Public License (Version 2) as 10 * published by the Free Software Foundation. 11 * 12 * This program is distributed in the hope it will be useful, but WITHOUT 13 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 14 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 15 * for more details. 16 * 17 * You should have received a copy of the GNU General Public License along 18 * with this program; if not, write to the Free Software Foundation, Inc., 19 * 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. 20 */ 21 22 #include <linux/compiler.h> 23 24 #include "ieee754dp.h" 25 26 int ieee754dp_class(union ieee754dp x) 27 { 28 COMPXDP; 29 EXPLODEXDP; 30 return xc; 31 } 32 33 static inline int ieee754dp_isnan(union ieee754dp x) 34 { 35 return ieee754_class_nan(ieee754dp_class(x)); 36 } 37 38 static inline int ieee754dp_issnan(union ieee754dp x) 39 { 40 assert(ieee754dp_isnan(x)); 41 return (DPMANT(x) & DP_MBIT(DP_FBITS - 1)) == DP_MBIT(DP_FBITS - 1); 42 } 43 44 45 /* 46 * Raise the Invalid Operation IEEE 754 exception 47 * and convert the signaling NaN supplied to a quiet NaN. 48 */ 49 union ieee754dp __cold ieee754dp_nanxcpt(union ieee754dp r) 50 { 51 assert(ieee754dp_issnan(r)); 52 53 ieee754_setcx(IEEE754_INVALID_OPERATION); 54 return ieee754dp_indef(); 55 } 56 57 static u64 ieee754dp_get_rounding(int sn, u64 xm) 58 { 59 /* inexact must round of 3 bits 60 */ 61 if (xm & (DP_MBIT(3) - 1)) { 62 switch (ieee754_csr.rm) { 63 case FPU_CSR_RZ: 64 break; 65 case FPU_CSR_RN: 66 xm += 0x3 + ((xm >> 3) & 1); 67 /* xm += (xm&0x8)?0x4:0x3 */ 68 break; 69 case FPU_CSR_RU: /* toward +Infinity */ 70 if (!sn) /* ?? */ 71 xm += 0x8; 72 break; 73 case FPU_CSR_RD: /* toward -Infinity */ 74 if (sn) /* ?? */ 75 xm += 0x8; 76 break; 77 } 78 } 79 return xm; 80 } 81 82 83 /* generate a normal/denormal number with over,under handling 84 * sn is sign 85 * xe is an unbiased exponent 86 * xm is 3bit extended precision value. 87 */ 88 union ieee754dp ieee754dp_format(int sn, int xe, u64 xm) 89 { 90 assert(xm); /* we don't gen exact zeros (probably should) */ 91 92 assert((xm >> (DP_FBITS + 1 + 3)) == 0); /* no execess */ 93 assert(xm & (DP_HIDDEN_BIT << 3)); 94 95 if (xe < DP_EMIN) { 96 /* strip lower bits */ 97 int es = DP_EMIN - xe; 98 99 if (ieee754_csr.nod) { 100 ieee754_setcx(IEEE754_UNDERFLOW); 101 ieee754_setcx(IEEE754_INEXACT); 102 103 switch(ieee754_csr.rm) { 104 case FPU_CSR_RN: 105 case FPU_CSR_RZ: 106 return ieee754dp_zero(sn); 107 case FPU_CSR_RU: /* toward +Infinity */ 108 if (sn == 0) 109 return ieee754dp_min(0); 110 else 111 return ieee754dp_zero(1); 112 case FPU_CSR_RD: /* toward -Infinity */ 113 if (sn == 0) 114 return ieee754dp_zero(0); 115 else 116 return ieee754dp_min(1); 117 } 118 } 119 120 if (xe == DP_EMIN - 1 && 121 ieee754dp_get_rounding(sn, xm) >> (DP_FBITS + 1 + 3)) 122 { 123 /* Not tiny after rounding */ 124 ieee754_setcx(IEEE754_INEXACT); 125 xm = ieee754dp_get_rounding(sn, xm); 126 xm >>= 1; 127 /* Clear grs bits */ 128 xm &= ~(DP_MBIT(3) - 1); 129 xe++; 130 } 131 else { 132 /* sticky right shift es bits 133 */ 134 xm = XDPSRS(xm, es); 135 xe += es; 136 assert((xm & (DP_HIDDEN_BIT << 3)) == 0); 137 assert(xe == DP_EMIN); 138 } 139 } 140 if (xm & (DP_MBIT(3) - 1)) { 141 ieee754_setcx(IEEE754_INEXACT); 142 if ((xm & (DP_HIDDEN_BIT << 3)) == 0) { 143 ieee754_setcx(IEEE754_UNDERFLOW); 144 } 145 146 /* inexact must round of 3 bits 147 */ 148 xm = ieee754dp_get_rounding(sn, xm); 149 /* adjust exponent for rounding add overflowing 150 */ 151 if (xm >> (DP_FBITS + 3 + 1)) { 152 /* add causes mantissa overflow */ 153 xm >>= 1; 154 xe++; 155 } 156 } 157 /* strip grs bits */ 158 xm >>= 3; 159 160 assert((xm >> (DP_FBITS + 1)) == 0); /* no execess */ 161 assert(xe >= DP_EMIN); 162 163 if (xe > DP_EMAX) { 164 ieee754_setcx(IEEE754_OVERFLOW); 165 ieee754_setcx(IEEE754_INEXACT); 166 /* -O can be table indexed by (rm,sn) */ 167 switch (ieee754_csr.rm) { 168 case FPU_CSR_RN: 169 return ieee754dp_inf(sn); 170 case FPU_CSR_RZ: 171 return ieee754dp_max(sn); 172 case FPU_CSR_RU: /* toward +Infinity */ 173 if (sn == 0) 174 return ieee754dp_inf(0); 175 else 176 return ieee754dp_max(1); 177 case FPU_CSR_RD: /* toward -Infinity */ 178 if (sn == 0) 179 return ieee754dp_max(0); 180 else 181 return ieee754dp_inf(1); 182 } 183 } 184 /* gen norm/denorm/zero */ 185 186 if ((xm & DP_HIDDEN_BIT) == 0) { 187 /* we underflow (tiny/zero) */ 188 assert(xe == DP_EMIN); 189 if (ieee754_csr.mx & IEEE754_UNDERFLOW) 190 ieee754_setcx(IEEE754_UNDERFLOW); 191 return builddp(sn, DP_EMIN - 1 + DP_EBIAS, xm); 192 } else { 193 assert((xm >> (DP_FBITS + 1)) == 0); /* no execess */ 194 assert(xm & DP_HIDDEN_BIT); 195 196 return builddp(sn, xe + DP_EBIAS, xm & ~DP_HIDDEN_BIT); 197 } 198 } 199