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