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
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 }