fpu_add.c revision 1.9 1 1.9 isaki /* $NetBSD: fpu_add.c,v 1.9 2013/03/26 11:30:20 isaki Exp $ */
2 1.1 briggs
3 1.1 briggs /*
4 1.1 briggs * Copyright (c) 1992, 1993
5 1.1 briggs * The Regents of the University of California. All rights reserved.
6 1.1 briggs *
7 1.1 briggs * This software was developed by the Computer Systems Engineering group
8 1.1 briggs * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
9 1.1 briggs * contributed to Berkeley.
10 1.1 briggs *
11 1.1 briggs * All advertising materials mentioning features or use of this software
12 1.1 briggs * must display the following acknowledgement:
13 1.1 briggs * This product includes software developed by the University of
14 1.1 briggs * California, Lawrence Berkeley Laboratory.
15 1.1 briggs *
16 1.1 briggs * Redistribution and use in source and binary forms, with or without
17 1.1 briggs * modification, are permitted provided that the following conditions
18 1.1 briggs * are met:
19 1.1 briggs * 1. Redistributions of source code must retain the above copyright
20 1.1 briggs * notice, this list of conditions and the following disclaimer.
21 1.1 briggs * 2. Redistributions in binary form must reproduce the above copyright
22 1.1 briggs * notice, this list of conditions and the following disclaimer in the
23 1.1 briggs * documentation and/or other materials provided with the distribution.
24 1.5 agc * 3. Neither the name of the University nor the names of its contributors
25 1.1 briggs * may be used to endorse or promote products derived from this software
26 1.1 briggs * without specific prior written permission.
27 1.1 briggs *
28 1.1 briggs * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
29 1.1 briggs * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
30 1.1 briggs * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
31 1.1 briggs * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
32 1.1 briggs * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
33 1.1 briggs * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
34 1.1 briggs * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
35 1.1 briggs * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
36 1.1 briggs * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
37 1.1 briggs * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
38 1.1 briggs * SUCH DAMAGE.
39 1.1 briggs *
40 1.1 briggs * @(#)fpu_add.c 8.1 (Berkeley) 6/11/93
41 1.1 briggs */
42 1.1 briggs
43 1.1 briggs /*
44 1.1 briggs * Perform an FPU add (return x + y).
45 1.1 briggs *
46 1.1 briggs * To subtract, negate y and call add.
47 1.1 briggs */
48 1.4 lukem
49 1.4 lukem #include <sys/cdefs.h>
50 1.9 isaki __KERNEL_RCSID(0, "$NetBSD: fpu_add.c,v 1.9 2013/03/26 11:30:20 isaki Exp $");
51 1.1 briggs
52 1.1 briggs #include <sys/types.h>
53 1.2 briggs #include <sys/systm.h>
54 1.1 briggs
55 1.1 briggs #include <machine/reg.h>
56 1.1 briggs
57 1.1 briggs #include "fpu_arith.h"
58 1.1 briggs #include "fpu_emulate.h"
59 1.1 briggs
60 1.1 briggs struct fpn *
61 1.8 isaki fpu_add(struct fpemu *fe)
62 1.1 briggs {
63 1.8 isaki struct fpn *x = &fe->fe_f1, *y = &fe->fe_f2, *r;
64 1.9 isaki uint32_t r0, r1, r2;
65 1.8 isaki int rd;
66 1.1 briggs
67 1.1 briggs /*
68 1.1 briggs * Put the `heavier' operand on the right (see fpu_emu.h).
69 1.1 briggs * Then we will have one of the following cases, taken in the
70 1.1 briggs * following order:
71 1.1 briggs *
72 1.1 briggs * - y = NaN. Implied: if only one is a signalling NaN, y is.
73 1.1 briggs * The result is y.
74 1.1 briggs * - y = Inf. Implied: x != NaN (is 0, number, or Inf: the NaN
75 1.1 briggs * case was taken care of earlier).
76 1.1 briggs * If x = -y, the result is NaN. Otherwise the result
77 1.1 briggs * is y (an Inf of whichever sign).
78 1.1 briggs * - y is 0. Implied: x = 0.
79 1.1 briggs * If x and y differ in sign (one positive, one negative),
80 1.1 briggs * the result is +0 except when rounding to -Inf. If same:
81 1.1 briggs * +0 + +0 = +0; -0 + -0 = -0.
82 1.1 briggs * - x is 0. Implied: y != 0.
83 1.1 briggs * Result is y.
84 1.1 briggs * - other. Implied: both x and y are numbers.
85 1.1 briggs * Do addition a la Hennessey & Patterson.
86 1.1 briggs */
87 1.1 briggs ORDER(x, y);
88 1.1 briggs if (ISNAN(y))
89 1.1 briggs return (y);
90 1.1 briggs if (ISINF(y)) {
91 1.1 briggs if (ISINF(x) && x->fp_sign != y->fp_sign)
92 1.1 briggs return (fpu_newnan(fe));
93 1.1 briggs return (y);
94 1.1 briggs }
95 1.1 briggs rd = (fe->fe_fpcr & FPCR_ROUND);
96 1.1 briggs if (ISZERO(y)) {
97 1.1 briggs if (rd != FPCR_MINF) /* only -0 + -0 gives -0 */
98 1.1 briggs y->fp_sign &= x->fp_sign;
99 1.1 briggs else /* any -0 operand gives -0 */
100 1.1 briggs y->fp_sign |= x->fp_sign;
101 1.1 briggs return (y);
102 1.1 briggs }
103 1.1 briggs if (ISZERO(x))
104 1.1 briggs return (y);
105 1.1 briggs /*
106 1.1 briggs * We really have two numbers to add, although their signs may
107 1.1 briggs * differ. Make the exponents match, by shifting the smaller
108 1.1 briggs * number right (e.g., 1.011 => 0.1011) and increasing its
109 1.1 briggs * exponent (2^3 => 2^4). Note that we do not alter the exponents
110 1.1 briggs * of x and y here.
111 1.1 briggs */
112 1.1 briggs r = &fe->fe_f3;
113 1.1 briggs r->fp_class = FPC_NUM;
114 1.1 briggs if (x->fp_exp == y->fp_exp) {
115 1.1 briggs r->fp_exp = x->fp_exp;
116 1.1 briggs r->fp_sticky = 0;
117 1.1 briggs } else {
118 1.1 briggs if (x->fp_exp < y->fp_exp) {
119 1.1 briggs /*
120 1.1 briggs * Try to avoid subtract case iii (see below).
121 1.1 briggs * This also guarantees that x->fp_sticky = 0.
122 1.1 briggs */
123 1.1 briggs SWAP(x, y);
124 1.1 briggs }
125 1.1 briggs /* now x->fp_exp > y->fp_exp */
126 1.1 briggs r->fp_exp = x->fp_exp;
127 1.1 briggs r->fp_sticky = fpu_shr(y, x->fp_exp - y->fp_exp);
128 1.1 briggs }
129 1.1 briggs r->fp_sign = x->fp_sign;
130 1.1 briggs if (x->fp_sign == y->fp_sign) {
131 1.1 briggs FPU_DECL_CARRY
132 1.1 briggs
133 1.1 briggs /*
134 1.1 briggs * The signs match, so we simply add the numbers. The result
135 1.1 briggs * may be `supernormal' (as big as 1.111...1 + 1.111...1, or
136 1.1 briggs * 11.111...0). If so, a single bit shift-right will fix it
137 1.1 briggs * (but remember to adjust the exponent).
138 1.1 briggs */
139 1.1 briggs /* r->fp_mant = x->fp_mant + y->fp_mant */
140 1.3 briggs FPU_ADDS(r->fp_mant[2], x->fp_mant[2], y->fp_mant[2]);
141 1.1 briggs FPU_ADDCS(r->fp_mant[1], x->fp_mant[1], y->fp_mant[1]);
142 1.1 briggs FPU_ADDC(r0, x->fp_mant[0], y->fp_mant[0]);
143 1.1 briggs if ((r->fp_mant[0] = r0) >= FP_2) {
144 1.1 briggs (void) fpu_shr(r, 1);
145 1.1 briggs r->fp_exp++;
146 1.1 briggs }
147 1.1 briggs } else {
148 1.1 briggs FPU_DECL_CARRY
149 1.1 briggs
150 1.1 briggs /*
151 1.1 briggs * The signs differ, so things are rather more difficult.
152 1.1 briggs * H&P would have us negate the negative operand and add;
153 1.1 briggs * this is the same as subtracting the negative operand.
154 1.1 briggs * This is quite a headache. Instead, we will subtract
155 1.1 briggs * y from x, regardless of whether y itself is the negative
156 1.1 briggs * operand. When this is done one of three conditions will
157 1.1 briggs * hold, depending on the magnitudes of x and y:
158 1.1 briggs * case i) |x| > |y|. The result is just x - y,
159 1.1 briggs * with x's sign, but it may need to be normalized.
160 1.1 briggs * case ii) |x| = |y|. The result is 0 (maybe -0)
161 1.1 briggs * so must be fixed up.
162 1.1 briggs * case iii) |x| < |y|. We goofed; the result should
163 1.1 briggs * be (y - x), with the same sign as y.
164 1.1 briggs * We could compare |x| and |y| here and avoid case iii,
165 1.1 briggs * but that would take just as much work as the subtract.
166 1.1 briggs * We can tell case iii has occurred by an overflow.
167 1.1 briggs *
168 1.1 briggs * N.B.: since x->fp_exp >= y->fp_exp, x->fp_sticky = 0.
169 1.1 briggs */
170 1.1 briggs /* r->fp_mant = x->fp_mant - y->fp_mant */
171 1.1 briggs FPU_SET_CARRY(y->fp_sticky);
172 1.1 briggs FPU_SUBCS(r2, x->fp_mant[2], y->fp_mant[2]);
173 1.1 briggs FPU_SUBCS(r1, x->fp_mant[1], y->fp_mant[1]);
174 1.1 briggs FPU_SUBC(r0, x->fp_mant[0], y->fp_mant[0]);
175 1.1 briggs if (r0 < FP_2) {
176 1.1 briggs /* cases i and ii */
177 1.3 briggs if ((r0 | r1 | r2) == 0) {
178 1.1 briggs /* case ii */
179 1.1 briggs r->fp_class = FPC_ZERO;
180 1.1 briggs r->fp_sign = (rd == FPCR_MINF);
181 1.1 briggs return (r);
182 1.1 briggs }
183 1.1 briggs } else {
184 1.1 briggs /*
185 1.1 briggs * Oops, case iii. This can only occur when the
186 1.1 briggs * exponents were equal, in which case neither
187 1.1 briggs * x nor y have sticky bits set. Flip the sign
188 1.1 briggs * (to y's sign) and negate the result to get y - x.
189 1.1 briggs */
190 1.1 briggs #ifdef DIAGNOSTIC
191 1.1 briggs if (x->fp_exp != y->fp_exp || r->fp_sticky)
192 1.1 briggs panic("fpu_add");
193 1.1 briggs #endif
194 1.1 briggs r->fp_sign = y->fp_sign;
195 1.3 briggs FPU_SUBS(r2, 0, r2);
196 1.1 briggs FPU_SUBCS(r1, 0, r1);
197 1.1 briggs FPU_SUBC(r0, 0, r0);
198 1.1 briggs }
199 1.1 briggs r->fp_mant[2] = r2;
200 1.1 briggs r->fp_mant[1] = r1;
201 1.1 briggs r->fp_mant[0] = r0;
202 1.1 briggs if (r0 < FP_1)
203 1.1 briggs fpu_norm(r);
204 1.1 briggs }
205 1.1 briggs return (r);
206 1.1 briggs }
207