fpu_rem.c revision 1.8 1 1.8 dsl /* $NetBSD: fpu_rem.c,v 1.8 2009/03/14 15:36:09 dsl Exp $ */
2 1.1 briggs
3 1.1 briggs /*
4 1.1 briggs * Copyright (c) 1995 Ken Nakata
5 1.1 briggs * All rights reserved.
6 1.1 briggs *
7 1.1 briggs * Redistribution and use in source and binary forms, with or without
8 1.1 briggs * modification, are permitted provided that the following conditions
9 1.1 briggs * are met:
10 1.1 briggs * 1. Redistributions of source code must retain the above copyright
11 1.1 briggs * notice, this list of conditions and the following disclaimer.
12 1.1 briggs * 2. Redistributions in binary form must reproduce the above copyright
13 1.1 briggs * notice, this list of conditions and the following disclaimer in the
14 1.1 briggs * documentation and/or other materials provided with the distribution.
15 1.1 briggs * 3. Neither the name of the author nor the names of its contributors
16 1.1 briggs * may be used to endorse or promote products derived from this software
17 1.1 briggs * without specific prior written permission.
18 1.1 briggs *
19 1.1 briggs * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
20 1.1 briggs * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21 1.1 briggs * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22 1.1 briggs * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
23 1.1 briggs * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
24 1.1 briggs * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
25 1.1 briggs * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
26 1.1 briggs * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
27 1.1 briggs * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
28 1.1 briggs * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
29 1.1 briggs * SUCH DAMAGE.
30 1.1 briggs *
31 1.1 briggs * @(#)fpu_rem.c 10/24/95
32 1.1 briggs */
33 1.5 lukem
34 1.5 lukem #include <sys/cdefs.h>
35 1.8 dsl __KERNEL_RCSID(0, "$NetBSD: fpu_rem.c,v 1.8 2009/03/14 15:36:09 dsl Exp $");
36 1.1 briggs
37 1.1 briggs #include <sys/types.h>
38 1.1 briggs #include <sys/signal.h>
39 1.1 briggs #include <machine/frame.h>
40 1.1 briggs
41 1.1 briggs #include "fpu_emulate.h"
42 1.1 briggs
43 1.1 briggs /*
44 1.1 briggs * ALGORITHM
45 1.1 briggs *
46 1.1 briggs * Step 1. Save and strip signs of X and Y: signX := sign(X),
47 1.1 briggs * signY := sign(Y), X := *X*, Y := *Y*,
48 1.1 briggs * signQ := signX EOR signY. Record whether MOD or REM
49 1.1 briggs * is requested.
50 1.1 briggs *
51 1.1 briggs * Step 2. Set L := expo(X)-expo(Y), k := 0, Q := 0.
52 1.1 briggs * If (L < 0) then
53 1.1 briggs * R := X, go to Step 4.
54 1.1 briggs * else
55 1.1 briggs * R := 2^(-L)X, j := L.
56 1.1 briggs * endif
57 1.1 briggs *
58 1.1 briggs * Step 3. Perform MOD(X,Y)
59 1.1 briggs * 3.1 If R = Y, go to Step 9.
60 1.1 briggs * 3.2 If R > Y, then { R := R - Y, Q := Q + 1}
61 1.1 briggs * 3.3 If j = 0, go to Step 4.
62 1.1 briggs * 3.4 k := k + 1, j := j - 1, Q := 2Q, R := 2R. Go to
63 1.1 briggs * Step 3.1.
64 1.1 briggs *
65 1.1 briggs * Step 4. At this point, R = X - QY = MOD(X,Y). Set
66 1.1 briggs * Last_Subtract := false (used in Step 7 below). If
67 1.1 briggs * MOD is requested, go to Step 6.
68 1.1 briggs *
69 1.1 briggs * Step 5. R = MOD(X,Y), but REM(X,Y) is requested.
70 1.1 briggs * 5.1 If R < Y/2, then R = MOD(X,Y) = REM(X,Y). Go to
71 1.1 briggs * Step 6.
72 1.1 briggs * 5.2 If R > Y/2, then { set Last_Subtract := true,
73 1.1 briggs * Q := Q + 1, Y := signY*Y }. Go to Step 6.
74 1.1 briggs * 5.3 This is the tricky case of R = Y/2. If Q is odd,
75 1.1 briggs * then { Q := Q + 1, signX := -signX }.
76 1.1 briggs *
77 1.1 briggs * Step 6. R := signX*R.
78 1.1 briggs *
79 1.1 briggs * Step 7. If Last_Subtract = true, R := R - Y.
80 1.1 briggs *
81 1.1 briggs * Step 8. Return signQ, last 7 bits of Q, and R as required.
82 1.1 briggs *
83 1.1 briggs * Step 9. At this point, R = 2^(-j)*X - Q Y = Y. Thus,
84 1.1 briggs * X = 2^(j)*(Q+1)Y. set Q := 2^(j)*(Q+1),
85 1.1 briggs * R := 0. Return signQ, last 7 bits of Q, and R.
86 1.1 briggs */
87 1.1 briggs
88 1.7 dsl static struct fpn * __fpu_modrem(struct fpemu *fe, int modrem);
89 1.1 briggs
90 1.1 briggs static struct fpn *
91 1.1 briggs __fpu_modrem(fe, modrem)
92 1.1 briggs struct fpemu *fe;
93 1.1 briggs int modrem;
94 1.1 briggs {
95 1.1 briggs static struct fpn X, Y;
96 1.1 briggs struct fpn *x, *y, *r;
97 1.1 briggs u_int signX, signY, signQ;
98 1.2 briggs int j, k, l, q;
99 1.1 briggs int Last_Subtract;
100 1.1 briggs
101 1.1 briggs CPYFPN(&X, &fe->fe_f1);
102 1.1 briggs CPYFPN(&Y, &fe->fe_f2);
103 1.1 briggs x = &X;
104 1.1 briggs y = &Y;
105 1.1 briggs r = &fe->fe_f2;
106 1.1 briggs
107 1.1 briggs /*
108 1.1 briggs * Step 1
109 1.1 briggs */
110 1.1 briggs signX = x->fp_sign;
111 1.1 briggs signY = y->fp_sign;
112 1.1 briggs signQ = (signX ^ signY);
113 1.1 briggs x->fp_sign = y->fp_sign = 0;
114 1.1 briggs
115 1.1 briggs /*
116 1.1 briggs * Step 2
117 1.1 briggs */
118 1.1 briggs l = x->fp_exp - y->fp_exp;
119 1.1 briggs k = 0;
120 1.1 briggs q = 0;
121 1.3 briggs if (l >= 0) {
122 1.1 briggs CPYFPN(r, x);
123 1.1 briggs r->fp_exp -= l;
124 1.1 briggs j = l;
125 1.1 briggs
126 1.1 briggs /*
127 1.1 briggs * Step 3
128 1.1 briggs */
129 1.1 briggs while (y->fp_exp != r->fp_exp || y->fp_mant[0] != r->fp_mant[0] ||
130 1.1 briggs y->fp_mant[1] != r->fp_mant[1] ||
131 1.4 briggs y->fp_mant[2] != r->fp_mant[2]) {
132 1.1 briggs
133 1.1 briggs /* Step 3.2 */
134 1.1 briggs if (y->fp_exp < r->fp_exp || y->fp_mant[0] < r->fp_mant[0] ||
135 1.1 briggs y->fp_mant[1] < r->fp_mant[1] ||
136 1.4 briggs y->fp_mant[2] < r->fp_mant[2]) {
137 1.1 briggs CPYFPN(&fe->fe_f1, r);
138 1.1 briggs CPYFPN(&fe->fe_f2, y);
139 1.1 briggs fe->fe_f2.fp_sign = 1;
140 1.1 briggs r = fpu_add(fe);
141 1.1 briggs q++;
142 1.1 briggs }
143 1.1 briggs
144 1.1 briggs /* Step 3.3 */
145 1.1 briggs if (j == 0)
146 1.1 briggs goto Step4;
147 1.1 briggs
148 1.1 briggs /* Step 3.4 */
149 1.1 briggs k++;
150 1.1 briggs j--;
151 1.1 briggs q += q;
152 1.1 briggs r->fp_exp++;
153 1.1 briggs }
154 1.3 briggs /* Step 9 */
155 1.3 briggs goto Step9;
156 1.1 briggs }
157 1.1 briggs Step4:
158 1.1 briggs Last_Subtract = 0;
159 1.1 briggs if (modrem == 0)
160 1.1 briggs goto Step6;
161 1.1 briggs
162 1.1 briggs /*
163 1.1 briggs * Step 5
164 1.1 briggs */
165 1.1 briggs /* Step 5.1 */
166 1.1 briggs if (r->fp_exp + 1 < y->fp_exp ||
167 1.2 briggs (r->fp_exp + 1 == y->fp_exp &&
168 1.2 briggs (r->fp_mant[0] < y->fp_mant[0] || r->fp_mant[1] < y->fp_mant[1] ||
169 1.4 briggs r->fp_mant[2] < y->fp_mant[2])))
170 1.1 briggs /* if r < y/2 */
171 1.1 briggs goto Step6;
172 1.1 briggs /* Step 5.2 */
173 1.1 briggs if (r->fp_exp + 1 != y->fp_exp ||
174 1.1 briggs r->fp_mant[0] != y->fp_mant[0] || r->fp_mant[1] != y->fp_mant[1] ||
175 1.4 briggs r->fp_mant[2] != y->fp_mant[2]) {
176 1.1 briggs /* if (!(r < y/2) && !(r == y/2)) */
177 1.1 briggs Last_Subtract = 1;
178 1.1 briggs q++;
179 1.1 briggs y->fp_sign = signY;
180 1.1 briggs } else {
181 1.1 briggs /* Step 5.3 */
182 1.1 briggs /* r == y/2 */
183 1.1 briggs if (q % 2) {
184 1.1 briggs q++;
185 1.1 briggs signX = !signX;
186 1.1 briggs }
187 1.1 briggs }
188 1.1 briggs
189 1.1 briggs Step6:
190 1.1 briggs r->fp_sign = signX;
191 1.1 briggs
192 1.1 briggs /*
193 1.1 briggs * Step 7
194 1.1 briggs */
195 1.1 briggs if (Last_Subtract) {
196 1.1 briggs CPYFPN(&fe->fe_f1, r);
197 1.1 briggs CPYFPN(&fe->fe_f2, y);
198 1.1 briggs fe->fe_f2.fp_sign = !y->fp_sign;
199 1.1 briggs r = fpu_add(fe);
200 1.1 briggs }
201 1.1 briggs /*
202 1.1 briggs * Step 8
203 1.1 briggs */
204 1.1 briggs q &= 0x7f;
205 1.1 briggs q |= (signQ << 7);
206 1.1 briggs fe->fe_fpframe->fpf_fpsr =
207 1.1 briggs fe->fe_fpsr =
208 1.1 briggs (fe->fe_fpsr & ~FPSR_QTT) | (q << 16);
209 1.1 briggs return r;
210 1.3 briggs
211 1.3 briggs Step9:
212 1.3 briggs fe->fe_f1.fp_class = FPC_ZERO;
213 1.3 briggs q++;
214 1.3 briggs q &= 0x7f;
215 1.3 briggs q |= (signQ << 7);
216 1.3 briggs fe->fe_fpframe->fpf_fpsr =
217 1.3 briggs fe->fe_fpsr =
218 1.3 briggs (fe->fe_fpsr & ~FPSR_QTT) | (q << 16);
219 1.3 briggs return &fe->fe_f1;
220 1.1 briggs }
221 1.1 briggs
222 1.1 briggs struct fpn *
223 1.8 dsl fpu_rem(struct fpemu *fe)
224 1.1 briggs {
225 1.1 briggs return __fpu_modrem(fe, 1);
226 1.1 briggs }
227 1.1 briggs
228 1.1 briggs struct fpn *
229 1.8 dsl fpu_mod(struct fpemu *fe)
230 1.1 briggs {
231 1.1 briggs return __fpu_modrem(fe, 0);
232 1.1 briggs }
233