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