PageRenderTime 47ms CodeModel.GetById 13ms app.highlight 27ms RepoModel.GetById 2ms app.codeStats 0ms

/arch/parisc/math-emu/dfrem.c

https://github.com/aicjofs/android_kernel_lge_v500_20d_f2fs
C | 297 lines | 149 code | 14 blank | 134 comment | 47 complexity | 3f748efc3c2abe1e256169e91e6a0d82 MD5 | raw file
  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
 50int
 51dbl_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}