1 // SPDX-License-Identifier: GPL-2.0-or-later 2 /* 3 * Linux/PA-RISC Project (http://www.parisc-linux.org/) 4 * 5 * Floating-point emulation code 6 * Copyright (C) 2001 Hewlett-Packard (Paul Bame) <bame@debian.org> 7 */ 8 /* 9 * BEGIN_DESC 10 * 11 * File: 12 * @(#) pa/spmath/sfmpy.c $Revision: 1.1 $ 13 * 14 * Purpose: 15 * Single Precision Floating-point Multiply 16 * 17 * External Interfaces: 18 * sgl_fmpy(srcptr1,srcptr2,dstptr,status) 19 * 20 * Internal Interfaces: 21 * 22 * Theory: 23 * <<please update with a overview of the operation of this file>> 24 * 25 * END_DESC 26 */ 27 28 29 #include "float.h" 30 #include "sgl_float.h" 31 32 /* 33 * Single Precision Floating-point Multiply 34 */ 35 36 int 37 sgl_fmpy( 38 sgl_floating_point *srcptr1, 39 sgl_floating_point *srcptr2, 40 sgl_floating_point *dstptr, 41 unsigned int *status) 42 { 43 register unsigned int opnd1, opnd2, opnd3, result; 44 register int dest_exponent, count; 45 register boolean inexact = FALSE, guardbit = FALSE, stickybit = FALSE; 46 boolean is_tiny; 47 48 opnd1 = *srcptr1; 49 opnd2 = *srcptr2; 50 /* 51 * set sign bit of result 52 */ 53 if (Sgl_sign(opnd1) ^ Sgl_sign(opnd2)) Sgl_setnegativezero(result); 54 else Sgl_setzero(result); 55 /* 56 * check first operand for NaN's or infinity 57 */ 58 if (Sgl_isinfinity_exponent(opnd1)) { 59 if (Sgl_iszero_mantissa(opnd1)) { 60 if (Sgl_isnotnan(opnd2)) { 61 if (Sgl_iszero_exponentmantissa(opnd2)) { 62 /* 63 * invalid since operands are infinity 64 * and zero 65 */ 66 if (Is_invalidtrap_enabled()) 67 return(INVALIDEXCEPTION); 68 Set_invalidflag(); 69 Sgl_makequietnan(result); 70 *dstptr = result; 71 return(NOEXCEPTION); 72 } 73 /* 74 * return infinity 75 */ 76 Sgl_setinfinity_exponentmantissa(result); 77 *dstptr = result; 78 return(NOEXCEPTION); 79 } 80 } 81 else { 82 /* 83 * is NaN; signaling or quiet? 84 */ 85 if (Sgl_isone_signaling(opnd1)) { 86 /* trap if INVALIDTRAP enabled */ 87 if (Is_invalidtrap_enabled()) 88 return(INVALIDEXCEPTION); 89 /* make NaN quiet */ 90 Set_invalidflag(); 91 Sgl_set_quiet(opnd1); 92 } 93 /* 94 * is second operand a signaling NaN? 95 */ 96 else if (Sgl_is_signalingnan(opnd2)) { 97 /* trap if INVALIDTRAP enabled */ 98 if (Is_invalidtrap_enabled()) 99 return(INVALIDEXCEPTION); 100 /* make NaN quiet */ 101 Set_invalidflag(); 102 Sgl_set_quiet(opnd2); 103 *dstptr = opnd2; 104 return(NOEXCEPTION); 105 } 106 /* 107 * return quiet NaN 108 */ 109 *dstptr = opnd1; 110 return(NOEXCEPTION); 111 } 112 } 113 /* 114 * check second operand for NaN's or infinity 115 */ 116 if (Sgl_isinfinity_exponent(opnd2)) { 117 if (Sgl_iszero_mantissa(opnd2)) { 118 if (Sgl_iszero_exponentmantissa(opnd1)) { 119 /* invalid since operands are zero & infinity */ 120 if (Is_invalidtrap_enabled()) 121 return(INVALIDEXCEPTION); 122 Set_invalidflag(); 123 Sgl_makequietnan(opnd2); 124 *dstptr = opnd2; 125 return(NOEXCEPTION); 126 } 127 /* 128 * return infinity 129 */ 130 Sgl_setinfinity_exponentmantissa(result); 131 *dstptr = result; 132 return(NOEXCEPTION); 133 } 134 /* 135 * is NaN; signaling or quiet? 136 */ 137 if (Sgl_isone_signaling(opnd2)) { 138 /* trap if INVALIDTRAP enabled */ 139 if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION); 140 141 /* make NaN quiet */ 142 Set_invalidflag(); 143 Sgl_set_quiet(opnd2); 144 } 145 /* 146 * return quiet NaN 147 */ 148 *dstptr = opnd2; 149 return(NOEXCEPTION); 150 } 151 /* 152 * Generate exponent 153 */ 154 dest_exponent = Sgl_exponent(opnd1) + Sgl_exponent(opnd2) - SGL_BIAS; 155 156 /* 157 * Generate mantissa 158 */ 159 if (Sgl_isnotzero_exponent(opnd1)) { 160 /* set hidden bit */ 161 Sgl_clear_signexponent_set_hidden(opnd1); 162 } 163 else { 164 /* check for zero */ 165 if (Sgl_iszero_mantissa(opnd1)) { 166 Sgl_setzero_exponentmantissa(result); 167 *dstptr = result; 168 return(NOEXCEPTION); 169 } 170 /* is denormalized, adjust exponent */ 171 Sgl_clear_signexponent(opnd1); 172 Sgl_leftshiftby1(opnd1); 173 Sgl_normalize(opnd1,dest_exponent); 174 } 175 /* opnd2 needs to have hidden bit set with msb in hidden bit */ 176 if (Sgl_isnotzero_exponent(opnd2)) { 177 Sgl_clear_signexponent_set_hidden(opnd2); 178 } 179 else { 180 /* check for zero */ 181 if (Sgl_iszero_mantissa(opnd2)) { 182 Sgl_setzero_exponentmantissa(result); 183 *dstptr = result; 184 return(NOEXCEPTION); 185 } 186 /* is denormalized; want to normalize */ 187 Sgl_clear_signexponent(opnd2); 188 Sgl_leftshiftby1(opnd2); 189 Sgl_normalize(opnd2,dest_exponent); 190 } 191 192 /* Multiply two source mantissas together */ 193 194 Sgl_leftshiftby4(opnd2); /* make room for guard bits */ 195 Sgl_setzero(opnd3); 196 /* 197 * Four bits at a time are inspected in each loop, and a 198 * simple shift and add multiply algorithm is used. 199 */ 200 for (count=1;count<SGL_P;count+=4) { 201 stickybit |= Slow4(opnd3); 202 Sgl_rightshiftby4(opnd3); 203 if (Sbit28(opnd1)) Sall(opnd3) += (Sall(opnd2) << 3); 204 if (Sbit29(opnd1)) Sall(opnd3) += (Sall(opnd2) << 2); 205 if (Sbit30(opnd1)) Sall(opnd3) += (Sall(opnd2) << 1); 206 if (Sbit31(opnd1)) Sall(opnd3) += Sall(opnd2); 207 Sgl_rightshiftby4(opnd1); 208 } 209 /* make sure result is left-justified */ 210 if (Sgl_iszero_sign(opnd3)) { 211 Sgl_leftshiftby1(opnd3); 212 } 213 else { 214 /* result mantissa >= 2. */ 215 dest_exponent++; 216 } 217 /* check for denormalized result */ 218 while (Sgl_iszero_sign(opnd3)) { 219 Sgl_leftshiftby1(opnd3); 220 dest_exponent--; 221 } 222 /* 223 * check for guard, sticky and inexact bits 224 */ 225 stickybit |= Sgl_all(opnd3) << (SGL_BITLENGTH - SGL_EXP_LENGTH + 1); 226 guardbit = Sbit24(opnd3); 227 inexact = guardbit | stickybit; 228 229 /* re-align mantissa */ 230 Sgl_rightshiftby8(opnd3); 231 232 /* 233 * round result 234 */ 235 if (inexact && (dest_exponent>0 || Is_underflowtrap_enabled())) { 236 Sgl_clear_signexponent(opnd3); 237 switch (Rounding_mode()) { 238 case ROUNDPLUS: 239 if (Sgl_iszero_sign(result)) 240 Sgl_increment(opnd3); 241 break; 242 case ROUNDMINUS: 243 if (Sgl_isone_sign(result)) 244 Sgl_increment(opnd3); 245 break; 246 case ROUNDNEAREST: 247 if (guardbit) { 248 if (stickybit || Sgl_isone_lowmantissa(opnd3)) 249 Sgl_increment(opnd3); 250 } 251 } 252 if (Sgl_isone_hidden(opnd3)) dest_exponent++; 253 } 254 Sgl_set_mantissa(result,opnd3); 255 256 /* 257 * Test for overflow 258 */ 259 if (dest_exponent >= SGL_INFINITY_EXPONENT) { 260 /* trap if OVERFLOWTRAP enabled */ 261 if (Is_overflowtrap_enabled()) { 262 /* 263 * Adjust bias of result 264 */ 265 Sgl_setwrapped_exponent(result,dest_exponent,ovfl); 266 *dstptr = result; 267 if (inexact) 268 if (Is_inexacttrap_enabled()) 269 return(OVERFLOWEXCEPTION | INEXACTEXCEPTION); 270 else Set_inexactflag(); 271 return(OVERFLOWEXCEPTION); 272 } 273 inexact = TRUE; 274 Set_overflowflag(); 275 /* set result to infinity or largest number */ 276 Sgl_setoverflow(result); 277 } 278 /* 279 * Test for underflow 280 */ 281 else if (dest_exponent <= 0) { 282 /* trap if UNDERFLOWTRAP enabled */ 283 if (Is_underflowtrap_enabled()) { 284 /* 285 * Adjust bias of result 286 */ 287 Sgl_setwrapped_exponent(result,dest_exponent,unfl); 288 *dstptr = result; 289 if (inexact) 290 if (Is_inexacttrap_enabled()) 291 return(UNDERFLOWEXCEPTION | INEXACTEXCEPTION); 292 else Set_inexactflag(); 293 return(UNDERFLOWEXCEPTION); 294 } 295 296 /* Determine if should set underflow flag */ 297 is_tiny = TRUE; 298 if (dest_exponent == 0 && inexact) { 299 switch (Rounding_mode()) { 300 case ROUNDPLUS: 301 if (Sgl_iszero_sign(result)) { 302 Sgl_increment(opnd3); 303 if (Sgl_isone_hiddenoverflow(opnd3)) 304 is_tiny = FALSE; 305 Sgl_decrement(opnd3); 306 } 307 break; 308 case ROUNDMINUS: 309 if (Sgl_isone_sign(result)) { 310 Sgl_increment(opnd3); 311 if (Sgl_isone_hiddenoverflow(opnd3)) 312 is_tiny = FALSE; 313 Sgl_decrement(opnd3); 314 } 315 break; 316 case ROUNDNEAREST: 317 if (guardbit && (stickybit || 318 Sgl_isone_lowmantissa(opnd3))) { 319 Sgl_increment(opnd3); 320 if (Sgl_isone_hiddenoverflow(opnd3)) 321 is_tiny = FALSE; 322 Sgl_decrement(opnd3); 323 } 324 break; 325 } 326 } 327 328 /* 329 * denormalize result or set to signed zero 330 */ 331 stickybit = inexact; 332 Sgl_denormalize(opnd3,dest_exponent,guardbit,stickybit,inexact); 333 334 /* return zero or smallest number */ 335 if (inexact) { 336 switch (Rounding_mode()) { 337 case ROUNDPLUS: 338 if (Sgl_iszero_sign(result)) { 339 Sgl_increment(opnd3); 340 } 341 break; 342 case ROUNDMINUS: 343 if (Sgl_isone_sign(result)) { 344 Sgl_increment(opnd3); 345 } 346 break; 347 case ROUNDNEAREST: 348 if (guardbit && (stickybit || 349 Sgl_isone_lowmantissa(opnd3))) { 350 Sgl_increment(opnd3); 351 } 352 break; 353 } 354 if (is_tiny) Set_underflowflag(); 355 } 356 Sgl_set_exponentmantissa(result,opnd3); 357 } 358 else Sgl_set_exponent(result,dest_exponent); 359 *dstptr = result; 360 361 /* check for inexact */ 362 if (inexact) { 363 if (Is_inexacttrap_enabled()) return(INEXACTEXCEPTION); 364 else Set_inexactflag(); 365 } 366 return(NOEXCEPTION); 367 } 368