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