1 /*
2 * CDDL HEADER START
3 *
4 * The contents of this file are subject to the terms of the
5 * Common Development and Distribution License, Version 1.0 only
6 * (the "License"). You may not use this file except in compliance
7 * with the License.
8 *
9 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
10 * or http://www.opensolaris.org/os/licensing.
11 * See the License for the specific language governing permissions
12 * and limitations under the License.
13 *
14 * When distributing Covered Code, include this CDDL HEADER in each
15 * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
16 * If applicable, add the following below this CDDL HEADER, with the
17 * fields enclosed by brackets "[]" replaced with your own identifying
18 * information: Portions Copyright [yyyy] [name of copyright owner]
19 *
20 * CDDL HEADER END
21 */
22 #pragma ident "%Z%%M% %I% %E% SMI"
23
24 /*
25 * Copyright (c) 1988 by Sun Microsystems, Inc.
26 */
27
28 #include "_Qquad.h"
29 #include "_Qglobals.h"
30
31 PRIVATE void
true_add(px,py,pz)32 true_add(px, py, pz)
33 unpacked *px, *py, *pz;
34
35 {
36 unsigned c;
37 unpacked *pt;
38
39 if ((int) px->fpclass < (int) py->fpclass) { /* Reverse. */
40 pt = py;
41 py = px;
42 px = pt;
43 }
44 /* Now class(x) >= class(y). */
45 switch (px->fpclass) {
46 case fp_quiet: /* NaN + x -> NaN */
47 case fp_signaling: /* NaN + x -> NaN */
48 case fp_infinity: /* Inf + x -> Inf */
49 case fp_zero: /* 0 + 0 -> 0 */
50 *pz = *px;
51 return;
52 default:
53 if (py->fpclass == fp_zero) {
54 *pz = *px;
55 return;
56 }
57 }
58 /* Now z is normal or subnormal. */
59 /* Now y is normal or subnormal. */
60 if (px->exponent < py->exponent) { /* Reverse. */
61 pt = py;
62 py = px;
63 px = pt;
64 }
65 /* Now class(x) >= class(y). */
66 pz->fpclass = px->fpclass;
67 pz->sign = px->sign;
68 pz->exponent = px->exponent;
69 pz->rounded = pz->sticky = 0;
70
71 if (px->exponent != py->exponent) { /* pre-alignment required */
72 fpu_rightshift(py, pz->exponent - py->exponent);
73 pz->rounded = py->rounded;
74 pz->sticky = py->sticky;
75 }
76 c = 0;
77 c = fpu_add3wc(&(pz->significand[3]),px->significand[3],
78 py->significand[3],c);
79 c = fpu_add3wc(&(pz->significand[2]),px->significand[2],
80 py->significand[2],c);
81 c = fpu_add3wc(&(pz->significand[1]),px->significand[1],
82 py->significand[1],c);
83 c = fpu_add3wc(&(pz->significand[0]),px->significand[0],
84 py->significand[0],c);
85
86 /* Handle carry out of msb. */
87 if(pz->significand[0]>=0x20000) {
88 fpu_rightshift(pz, 1); /* Carried out bit. */
89 pz->exponent ++; /* Renormalize. */
90 }
91 return;
92 }
93
94 PRIVATE void
true_sub(px,py,pz)95 true_sub(px, py, pz)
96 unpacked *px, *py, *pz;
97
98 {
99 unsigned *z,g,s,r,c;
100 int n;
101 unpacked *pt;
102
103 if ((int) px->fpclass < (int) py->fpclass) { /* Reverse. */
104 pt = py;
105 py = px;
106 px = pt;
107 }
108 /* Now class(x) >= class(y). */
109 *pz = *px; /* Tentative difference: x. */
110 switch (pz->fpclass) {
111 case fp_quiet: /* NaN - x -> NaN */
112 case fp_signaling: /* NaN - x -> NaN */
113 return;
114 case fp_infinity: /* Inf - x -> Inf */
115 if (py->fpclass == fp_infinity) {
116 fpu_error_nan(pz); /* Inf - Inf -> NaN */
117 pz->fpclass = fp_quiet;
118 }
119 return;
120 case fp_zero: /* 0 - 0 -> 0 */
121 pz->sign = (fp_direction == fp_negative);
122 return;
123 default:
124 if (py->fpclass == fp_zero)
125 return;
126 }
127
128 /* x and y are both normal or subnormal. */
129
130 if (px->exponent < py->exponent) { /* Reverse. */
131 pt = py;
132 py = px;
133 px = pt;
134 }
135 /* Now exp(x) >= exp(y). */
136 pz->fpclass = px->fpclass;
137 pz->sign = px->sign;
138 pz->exponent = px->exponent;
139 pz->rounded = 0;
140 pz->sticky = 0;
141 z = pz->significand;
142
143 if (px->exponent == py->exponent) { /* no pre-alignment required */
144 c = 0;
145 c = fpu_sub3wc(&z[3],px->significand[3],py->significand[3],c);
146 c = fpu_sub3wc(&z[2],px->significand[2],py->significand[2],c);
147 c = fpu_sub3wc(&z[1],px->significand[1],py->significand[1],c);
148 c = fpu_sub3wc(&z[0],px->significand[0],py->significand[0],c);
149 if((z[0]|z[1]|z[2]|z[3])==0) { /* exact zero result */
150 pz->sign = (fp_direction == fp_negative);
151 pz->fpclass = fp_zero;
152 return;
153 }
154 if(z[0]>=0x20000) { /* sign reversal occurred */
155 pz->sign = py->sign;
156 c = 0;
157 c = fpu_neg2wc(&z[3],z[3],c);
158 c = fpu_neg2wc(&z[2],z[2],c);
159 c = fpu_neg2wc(&z[1],z[1],c);
160 c = fpu_neg2wc(&z[0],z[0],c);
161 }
162 fpu_normalize(pz);
163 return;
164 } else { /* pre-alignment required */
165 fpu_rightshift(py, pz->exponent - py->exponent - 1);
166 r = py->rounded; /* rounded bit */
167 s = py->sticky; /* sticky bit */
168 fpu_rightshift(py, 1);
169 g = py->rounded; /* guard bit */
170 if(s!=0) r = (r==0);
171 if((r|s)!=0) g = (g==0);/* guard and rounded bits of z */
172 c = ((g|r|s)!=0);
173 c = fpu_sub3wc(&z[3],px->significand[3],py->significand[3],c);
174 c = fpu_sub3wc(&z[2],px->significand[2],py->significand[2],c);
175 c = fpu_sub3wc(&z[1],px->significand[1],py->significand[1],c);
176 c = fpu_sub3wc(&z[0],px->significand[0],py->significand[0],c);
177
178 if(z[0]>=0x10000) { /* don't need post-shifted */
179 pz->sticky = s|r;
180 pz->rounded = g;
181 } else { /* post-shifted left 1 bit */
182 pz->sticky = s;
183 pz->rounded = r;
184 pz->significand[0] = (z[0]<<1)|((z[1]&0x80000000)>>31);
185 pz->significand[1] = (z[1]<<1)|((z[2]&0x80000000)>>31);
186 pz->significand[2] = (z[2]<<1)|((z[3]&0x80000000)>>31);
187 pz->significand[3] = (z[3]<<1)|g;
188 pz->exponent -= 1;
189 if(z[0]<0x10000) fpu_normalize(pz);
190 }
191 return;
192 }
193 }
194
195 void
_fp_add(px,py,pz)196 _fp_add(px, py, pz)
197 unpacked *px, *py, *pz;
198
199 {
200 if (px->sign == py->sign)
201 true_add(px, py, pz);
202 else
203 true_sub(px, py, pz);
204 }
205
206 void
_fp_sub(px,py,pz)207 _fp_sub(px, py, pz)
208 unpacked *px, *py, *pz;
209
210 {
211 py->sign = 1 - py->sign;
212 if (px->sign == py->sign)
213 true_add(px, py, pz);
214 else
215 true_sub(px, py, pz);
216 }
217