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
sgl_fmpy(sgl_floating_point * srcptr1,sgl_floating_point * srcptr2,sgl_floating_point * dstptr,unsigned int * status)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