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/sfrem.c		$Revision: 1.1 $
13  *
14  *  Purpose:
15  *	Single Precision Floating-point Remainder
16  *
17  *  External Interfaces:
18  *	sgl_frem(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 
30 #include "float.h"
31 #include "sgl_float.h"
32 
33 /*
34  *  Single Precision Floating-point Remainder
35  */
36 
37 int
sgl_frem(sgl_floating_point * srcptr1,sgl_floating_point * srcptr2,sgl_floating_point * dstptr,unsigned int * status)38 sgl_frem (sgl_floating_point * srcptr1, sgl_floating_point * srcptr2,
39 	  sgl_floating_point * dstptr, unsigned int *status)
40 {
41 	register unsigned int opnd1, opnd2, result;
42 	register int opnd1_exponent, opnd2_exponent, dest_exponent, stepcount;
43 	register boolean roundup = FALSE;
44 
45 	opnd1 = *srcptr1;
46 	opnd2 = *srcptr2;
47 	/*
48 	 * check first operand for NaN's or infinity
49 	 */
50 	if ((opnd1_exponent = Sgl_exponent(opnd1)) == SGL_INFINITY_EXPONENT) {
51 		if (Sgl_iszero_mantissa(opnd1)) {
52 			if (Sgl_isnotnan(opnd2)) {
53 				/* invalid since first operand is infinity */
54 				if (Is_invalidtrap_enabled())
55                                 	return(INVALIDEXCEPTION);
56                                 Set_invalidflag();
57                                 Sgl_makequietnan(result);
58 				*dstptr = result;
59 				return(NOEXCEPTION);
60 			}
61 		}
62 		else {
63                 	/*
64                  	 * is NaN; signaling or quiet?
65                  	 */
66                 	if (Sgl_isone_signaling(opnd1)) {
67                         	/* trap if INVALIDTRAP enabled */
68                         	if (Is_invalidtrap_enabled())
69                             		return(INVALIDEXCEPTION);
70                         	/* make NaN quiet */
71                         	Set_invalidflag();
72                         	Sgl_set_quiet(opnd1);
73                 	}
74 			/*
75 			 * is second operand a signaling NaN?
76 			 */
77 			else if (Sgl_is_signalingnan(opnd2)) {
78                         	/* trap if INVALIDTRAP enabled */
79                         	if (Is_invalidtrap_enabled())
80                             		return(INVALIDEXCEPTION);
81                         	/* make NaN quiet */
82                         	Set_invalidflag();
83                         	Sgl_set_quiet(opnd2);
84                 		*dstptr = opnd2;
85                 		return(NOEXCEPTION);
86 			}
87                 	/*
88                  	 * return quiet NaN
89                  	 */
90                 	*dstptr = opnd1;
91                 	return(NOEXCEPTION);
92 		}
93 	}
94 	/*
95 	 * check second operand for NaN's or infinity
96 	 */
97 	if ((opnd2_exponent = Sgl_exponent(opnd2)) == SGL_INFINITY_EXPONENT) {
98 		if (Sgl_iszero_mantissa(opnd2)) {
99 			/*
100 			 * return first operand
101 			 */
102                 	*dstptr = opnd1;
103 			return(NOEXCEPTION);
104 		}
105                 /*
106                  * is NaN; signaling or quiet?
107                  */
108                 if (Sgl_isone_signaling(opnd2)) {
109                         /* trap if INVALIDTRAP enabled */
110                         if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
111                         /* make NaN quiet */
112                         Set_invalidflag();
113                         Sgl_set_quiet(opnd2);
114                 }
115                 /*
116                  * return quiet NaN
117                  */
118                 *dstptr = opnd2;
119                 return(NOEXCEPTION);
120 	}
121 	/*
122 	 * check second operand for zero
123 	 */
124 	if (Sgl_iszero_exponentmantissa(opnd2)) {
125 		/* invalid since second operand is zero */
126 		if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
127                 Set_invalidflag();
128                 Sgl_makequietnan(result);
129 		*dstptr = result;
130 		return(NOEXCEPTION);
131 	}
132 
133 	/*
134 	 * get sign of result
135 	 */
136 	result = opnd1;
137 
138 	/*
139 	 * check for denormalized operands
140 	 */
141 	if (opnd1_exponent == 0) {
142 		/* check for zero */
143 		if (Sgl_iszero_mantissa(opnd1)) {
144 			*dstptr = opnd1;
145 			return(NOEXCEPTION);
146 		}
147 		/* normalize, then continue */
148 		opnd1_exponent = 1;
149 		Sgl_normalize(opnd1,opnd1_exponent);
150 	}
151 	else {
152 		Sgl_clear_signexponent_set_hidden(opnd1);
153 	}
154 	if (opnd2_exponent == 0) {
155 		/* normalize, then continue */
156 		opnd2_exponent = 1;
157 		Sgl_normalize(opnd2,opnd2_exponent);
158 	}
159 	else {
160 		Sgl_clear_signexponent_set_hidden(opnd2);
161 	}
162 
163 	/* find result exponent and divide step loop count */
164 	dest_exponent = opnd2_exponent - 1;
165 	stepcount = opnd1_exponent - opnd2_exponent;
166 
167 	/*
168 	 * check for opnd1/opnd2 < 1
169 	 */
170 	if (stepcount < 0) {
171 		/*
172 		 * check for opnd1/opnd2 > 1/2
173 		 *
174 		 * In this case n will round to 1, so
175 		 *    r = opnd1 - opnd2
176 		 */
177 		if (stepcount == -1 && Sgl_isgreaterthan(opnd1,opnd2)) {
178 			Sgl_all(result) = ~Sgl_all(result);   /* set sign */
179 			/* align opnd2 with opnd1 */
180 			Sgl_leftshiftby1(opnd2);
181 			Sgl_subtract(opnd2,opnd1,opnd2);
182 			/* now normalize */
183                 	while (Sgl_iszero_hidden(opnd2)) {
184                         	Sgl_leftshiftby1(opnd2);
185                         	dest_exponent--;
186 			}
187 			Sgl_set_exponentmantissa(result,opnd2);
188 			goto testforunderflow;
189 		}
190 		/*
191 		 * opnd1/opnd2 <= 1/2
192 		 *
193 		 * In this case n will round to zero, so
194 		 *    r = opnd1
195 		 */
196 		Sgl_set_exponentmantissa(result,opnd1);
197 		dest_exponent = opnd1_exponent;
198 		goto testforunderflow;
199 	}
200 
201 	/*
202 	 * Generate result
203 	 *
204 	 * Do iterative subtract until remainder is less than operand 2.
205 	 */
206 	while (stepcount-- > 0 && Sgl_all(opnd1)) {
207 		if (Sgl_isnotlessthan(opnd1,opnd2))
208 			Sgl_subtract(opnd1,opnd2,opnd1);
209 		Sgl_leftshiftby1(opnd1);
210 	}
211 	/*
212 	 * Do last subtract, then determine which way to round if remainder
213 	 * is exactly 1/2 of opnd2
214 	 */
215 	if (Sgl_isnotlessthan(opnd1,opnd2)) {
216 		Sgl_subtract(opnd1,opnd2,opnd1);
217 		roundup = TRUE;
218 	}
219 	if (stepcount > 0 || Sgl_iszero(opnd1)) {
220 		/* division is exact, remainder is zero */
221 		Sgl_setzero_exponentmantissa(result);
222 		*dstptr = result;
223 		return(NOEXCEPTION);
224 	}
225 
226 	/*
227 	 * Check for cases where opnd1/opnd2 < n
228 	 *
229 	 * In this case the result's sign will be opposite that of
230 	 * opnd1.  The mantissa also needs some correction.
231 	 */
232 	Sgl_leftshiftby1(opnd1);
233 	if (Sgl_isgreaterthan(opnd1,opnd2)) {
234 		Sgl_invert_sign(result);
235 		Sgl_subtract((opnd2<<1),opnd1,opnd1);
236 	}
237 	/* check for remainder being exactly 1/2 of opnd2 */
238 	else if (Sgl_isequal(opnd1,opnd2) && roundup) {
239 		Sgl_invert_sign(result);
240 	}
241 
242 	/* normalize result's mantissa */
243         while (Sgl_iszero_hidden(opnd1)) {
244                 dest_exponent--;
245                 Sgl_leftshiftby1(opnd1);
246         }
247 	Sgl_set_exponentmantissa(result,opnd1);
248 
249         /*
250          * Test for underflow
251          */
252     testforunderflow:
253 	if (dest_exponent <= 0) {
254                 /* trap if UNDERFLOWTRAP enabled */
255                 if (Is_underflowtrap_enabled()) {
256                         /*
257                          * Adjust bias of result
258                          */
259                         Sgl_setwrapped_exponent(result,dest_exponent,unfl);
260 			*dstptr = result;
261 			/* frem is always exact */
262 			return(UNDERFLOWEXCEPTION);
263                 }
264                 /*
265                  * denormalize result or set to signed zero
266                  */
267                 if (dest_exponent >= (1 - SGL_P)) {
268 			Sgl_rightshift_exponentmantissa(result,1-dest_exponent);
269                 }
270                 else {
271 			Sgl_setzero_exponentmantissa(result);
272 		}
273 	}
274 	else Sgl_set_exponent(result,dest_exponent);
275 	*dstptr = result;
276 	return(NOEXCEPTION);
277 }
278