ARM: 6979/1: mach-vt8500: add forgotten irq_data conversion
[linux-2.6.git] / arch / parisc / math-emu / dfrem.c
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/dfrem.c               $Revision: 1.1 $
26  *
27  *  Purpose:
28  *      Double Precision Floating-point Remainder
29  *
30  *  External Interfaces:
31  *      dbl_frem(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
43 #include "float.h"
44 #include "dbl_float.h"
45
46 /*
47  *  Double Precision Floating-point Remainder
48  */
49
50 int
51 dbl_frem (dbl_floating_point * srcptr1, dbl_floating_point * srcptr2,
52           dbl_floating_point * dstptr, unsigned int *status)
53 {
54         register unsigned int opnd1p1, opnd1p2, opnd2p1, opnd2p2;
55         register unsigned int resultp1, resultp2;
56         register int opnd1_exponent, opnd2_exponent, dest_exponent, stepcount;
57         register boolean roundup = FALSE;
58
59         Dbl_copyfromptr(srcptr1,opnd1p1,opnd1p2);
60         Dbl_copyfromptr(srcptr2,opnd2p1,opnd2p2);
61         /*
62          * check first operand for NaN's or infinity
63          */
64         if ((opnd1_exponent = Dbl_exponent(opnd1p1)) == DBL_INFINITY_EXPONENT) {
65                 if (Dbl_iszero_mantissa(opnd1p1,opnd1p2)) {
66                         if (Dbl_isnotnan(opnd2p1,opnd2p2)) {
67                                 /* invalid since first operand is infinity */
68                                 if (Is_invalidtrap_enabled()) 
69                                         return(INVALIDEXCEPTION);
70                                 Set_invalidflag();
71                                 Dbl_makequietnan(resultp1,resultp2);
72                                 Dbl_copytoptr(resultp1,resultp2,dstptr);
73                                 return(NOEXCEPTION);
74                         }
75                 }
76                 else {
77                         /*
78                          * is NaN; signaling or quiet?
79                          */
80                         if (Dbl_isone_signaling(opnd1p1)) {
81                                 /* trap if INVALIDTRAP enabled */
82                                 if (Is_invalidtrap_enabled()) 
83                                         return(INVALIDEXCEPTION);
84                                 /* make NaN quiet */
85                                 Set_invalidflag();
86                                 Dbl_set_quiet(opnd1p1);
87                         }
88                         /* 
89                          * is second operand a signaling NaN? 
90                          */
91                         else if (Dbl_is_signalingnan(opnd2p1)) {
92                                 /* trap if INVALIDTRAP enabled */
93                                 if (Is_invalidtrap_enabled()) 
94                                         return(INVALIDEXCEPTION);
95                                 /* make NaN quiet */
96                                 Set_invalidflag();
97                                 Dbl_set_quiet(opnd2p1);
98                                 Dbl_copytoptr(opnd2p1,opnd2p2,dstptr);
99                                 return(NOEXCEPTION);
100                         }
101                         /*
102                          * return quiet NaN
103                          */
104                         Dbl_copytoptr(opnd1p1,opnd1p2,dstptr);
105                         return(NOEXCEPTION);
106                 }
107         } 
108         /*
109          * check second operand for NaN's or infinity
110          */
111         if ((opnd2_exponent = Dbl_exponent(opnd2p1)) == DBL_INFINITY_EXPONENT) {
112                 if (Dbl_iszero_mantissa(opnd2p1,opnd2p2)) {
113                         /*
114                          * return first operand
115                          */
116                         Dbl_copytoptr(opnd1p1,opnd1p2,dstptr);
117                         return(NOEXCEPTION);
118                 }
119                 /*
120                  * is NaN; signaling or quiet?
121                  */
122                 if (Dbl_isone_signaling(opnd2p1)) {
123                         /* trap if INVALIDTRAP enabled */
124                         if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
125                         /* make NaN quiet */
126                         Set_invalidflag();
127                         Dbl_set_quiet(opnd2p1);
128                 }
129                 /*
130                  * return quiet NaN
131                  */
132                 Dbl_copytoptr(opnd2p1,opnd2p2,dstptr);
133                 return(NOEXCEPTION);
134         }
135         /*
136          * check second operand for zero
137          */
138         if (Dbl_iszero_exponentmantissa(opnd2p1,opnd2p2)) {
139                 /* invalid since second operand is zero */
140                 if (Is_invalidtrap_enabled()) return(INVALIDEXCEPTION);
141                 Set_invalidflag();
142                 Dbl_makequietnan(resultp1,resultp2);
143                 Dbl_copytoptr(resultp1,resultp2,dstptr);
144                 return(NOEXCEPTION);
145         }
146
147         /* 
148          * get sign of result
149          */
150         resultp1 = opnd1p1;  
151
152         /* 
153          * check for denormalized operands
154          */
155         if (opnd1_exponent == 0) {
156                 /* check for zero */
157                 if (Dbl_iszero_mantissa(opnd1p1,opnd1p2)) {
158                         Dbl_copytoptr(opnd1p1,opnd1p2,dstptr);
159                         return(NOEXCEPTION);
160                 }
161                 /* normalize, then continue */
162                 opnd1_exponent = 1;
163                 Dbl_normalize(opnd1p1,opnd1p2,opnd1_exponent);
164         }
165         else {
166                 Dbl_clear_signexponent_set_hidden(opnd1p1);
167         }
168         if (opnd2_exponent == 0) {
169                 /* normalize, then continue */
170                 opnd2_exponent = 1;
171                 Dbl_normalize(opnd2p1,opnd2p2,opnd2_exponent);
172         }
173         else {
174                 Dbl_clear_signexponent_set_hidden(opnd2p1);
175         }
176
177         /* find result exponent and divide step loop count */
178         dest_exponent = opnd2_exponent - 1;
179         stepcount = opnd1_exponent - opnd2_exponent;
180
181         /*
182          * check for opnd1/opnd2 < 1
183          */
184         if (stepcount < 0) {
185                 /*
186                  * check for opnd1/opnd2 > 1/2
187                  *
188                  * In this case n will round to 1, so 
189                  *    r = opnd1 - opnd2 
190                  */
191                 if (stepcount == -1 && 
192                     Dbl_isgreaterthan(opnd1p1,opnd1p2,opnd2p1,opnd2p2)) {
193                         /* set sign */
194                         Dbl_allp1(resultp1) = ~Dbl_allp1(resultp1);
195                         /* align opnd2 with opnd1 */
196                         Dbl_leftshiftby1(opnd2p1,opnd2p2); 
197                         Dbl_subtract(opnd2p1,opnd2p2,opnd1p1,opnd1p2,
198                          opnd2p1,opnd2p2);
199                         /* now normalize */
200                         while (Dbl_iszero_hidden(opnd2p1)) {
201                                 Dbl_leftshiftby1(opnd2p1,opnd2p2);
202                                 dest_exponent--;
203                         }
204                         Dbl_set_exponentmantissa(resultp1,resultp2,opnd2p1,opnd2p2);
205                         goto testforunderflow;
206                 }
207                 /*
208                  * opnd1/opnd2 <= 1/2
209                  *
210                  * In this case n will round to zero, so 
211                  *    r = opnd1
212                  */
213                 Dbl_set_exponentmantissa(resultp1,resultp2,opnd1p1,opnd1p2);
214                 dest_exponent = opnd1_exponent;
215                 goto testforunderflow;
216         }
217
218         /*
219          * Generate result
220          *
221          * Do iterative subtract until remainder is less than operand 2.
222          */
223         while (stepcount-- > 0 && (Dbl_allp1(opnd1p1) || Dbl_allp2(opnd1p2))) {
224                 if (Dbl_isnotlessthan(opnd1p1,opnd1p2,opnd2p1,opnd2p2)) {
225                         Dbl_subtract(opnd1p1,opnd1p2,opnd2p1,opnd2p2,opnd1p1,opnd1p2);
226                 }
227                 Dbl_leftshiftby1(opnd1p1,opnd1p2);
228         }
229         /*
230          * Do last subtract, then determine which way to round if remainder 
231          * is exactly 1/2 of opnd2 
232          */
233         if (Dbl_isnotlessthan(opnd1p1,opnd1p2,opnd2p1,opnd2p2)) {
234                 Dbl_subtract(opnd1p1,opnd1p2,opnd2p1,opnd2p2,opnd1p1,opnd1p2);
235                 roundup = TRUE;
236         }
237         if (stepcount > 0 || Dbl_iszero(opnd1p1,opnd1p2)) {
238                 /* division is exact, remainder is zero */
239                 Dbl_setzero_exponentmantissa(resultp1,resultp2);
240                 Dbl_copytoptr(resultp1,resultp2,dstptr);
241                 return(NOEXCEPTION);
242         }
243
244         /* 
245          * Check for cases where opnd1/opnd2 < n 
246          *
247          * In this case the result's sign will be opposite that of
248          * opnd1.  The mantissa also needs some correction.
249          */
250         Dbl_leftshiftby1(opnd1p1,opnd1p2);
251         if (Dbl_isgreaterthan(opnd1p1,opnd1p2,opnd2p1,opnd2p2)) {
252                 Dbl_invert_sign(resultp1);
253                 Dbl_leftshiftby1(opnd2p1,opnd2p2);
254                 Dbl_subtract(opnd2p1,opnd2p2,opnd1p1,opnd1p2,opnd1p1,opnd1p2);
255         }
256         /* check for remainder being exactly 1/2 of opnd2 */
257         else if (Dbl_isequal(opnd1p1,opnd1p2,opnd2p1,opnd2p2) && roundup) { 
258                 Dbl_invert_sign(resultp1);
259         }
260
261         /* normalize result's mantissa */
262         while (Dbl_iszero_hidden(opnd1p1)) {
263                 dest_exponent--;
264                 Dbl_leftshiftby1(opnd1p1,opnd1p2);
265         }
266         Dbl_set_exponentmantissa(resultp1,resultp2,opnd1p1,opnd1p2);
267
268         /* 
269          * Test for underflow
270          */
271     testforunderflow:
272         if (dest_exponent <= 0) {
273                 /* trap if UNDERFLOWTRAP enabled */
274                 if (Is_underflowtrap_enabled()) {
275                         /*
276                          * Adjust bias of result
277                          */
278                         Dbl_setwrapped_exponent(resultp1,dest_exponent,unfl);
279                         /* frem is always exact */
280                         Dbl_copytoptr(resultp1,resultp2,dstptr);
281                         return(UNDERFLOWEXCEPTION);
282                 }
283                 /*
284                  * denormalize result or set to signed zero
285                  */
286                 if (dest_exponent >= (1 - DBL_P)) {
287                         Dbl_rightshift_exponentmantissa(resultp1,resultp2,
288                          1-dest_exponent);
289                 }
290                 else {
291                         Dbl_setzero_exponentmantissa(resultp1,resultp2);
292                 }
293         }
294         else Dbl_set_exponent(resultp1,dest_exponent);
295         Dbl_copytoptr(resultp1,resultp2,dstptr);
296         return(NOEXCEPTION);
297 }