]>
Commit | Line | Data |
---|---|---|
252b5132 | 1 | /* flonum_mult.c - multiply two flonums |
2da5c037 | 2 | Copyright 1987, 1990, 1991, 1992, 1995, 2000, 2002, 2003 |
f0e652b4 | 3 | Free Software Foundation, Inc. |
252b5132 RH |
4 | |
5 | This file is part of Gas, the GNU Assembler. | |
6 | ||
7 | The GNU assembler is distributed in the hope that it will be | |
8 | useful, but WITHOUT ANY WARRANTY. No author or distributor | |
9 | accepts responsibility to anyone for the consequences of using it | |
10 | or for whether it serves any particular purpose or works at all, | |
11 | unless he says so in writing. Refer to the GNU Assembler General | |
12 | Public License for full details. | |
13 | ||
14 | Everyone is granted permission to copy, modify and redistribute | |
15 | the GNU Assembler, but only under the conditions described in the | |
16 | GNU Assembler General Public License. A copy of this license is | |
17 | supposed to have been given to you along with the GNU Assembler | |
18 | so you can know your rights and responsibilities. It should be | |
19 | in a file named COPYING. Among other things, the copyright | |
20 | notice and this notice must be preserved on all copies. */ | |
21 | ||
29589b0c | 22 | #include "ansidecl.h" |
252b5132 RH |
23 | #include "flonum.h" |
24 | ||
25 | /* plan for a . b => p(roduct) | |
f851444e | 26 | |
252b5132 RH |
27 | +-------+-------+-/ /-+-------+-------+ |
28 | | a | a | ... | a | a | | |
29 | | A | A-1 | | 1 | 0 | | |
30 | +-------+-------+-/ /-+-------+-------+ | |
f851444e | 31 | |
252b5132 RH |
32 | +-------+-------+-/ /-+-------+-------+ |
33 | | b | b | ... | b | b | | |
34 | | B | B-1 | | 1 | 0 | | |
35 | +-------+-------+-/ /-+-------+-------+ | |
f851444e | 36 | |
252b5132 RH |
37 | +-------+-------+-/ /-+-------+-/ /-+-------+-------+ |
38 | | p | p | ... | p | ... | p | p | | |
39 | | A+B+1| A+B | | N | | 1 | 0 | | |
40 | +-------+-------+-/ /-+-------+-/ /-+-------+-------+ | |
f851444e | 41 | |
252b5132 RH |
42 | /^\ |
43 | (carry) a .b ... | ... a .b a .b | |
44 | A B | 0 1 0 0 | |
45 | | | |
46 | ... | ... a .b | |
47 | | 1 0 | |
48 | | | |
49 | | ... | |
50 | | | |
51 | | | |
52 | | | |
53 | | ___ | |
54 | | \ | |
55 | +----- P = > a .b | |
56 | N /__ i j | |
f851444e | 57 | |
252b5132 | 58 | N = 0 ... A+B |
f851444e | 59 | |
252b5132 RH |
60 | for all i,j where i+j=N |
61 | [i,j integers > 0] | |
f851444e | 62 | |
252b5132 RH |
63 | a[], b[], p[] may not intersect. |
64 | Zero length factors signify 0 significant bits: treat as 0.0. | |
65 | 0.0 factors do the right thing. | |
66 | Zero length product OK. | |
f851444e | 67 | |
252b5132 RH |
68 | I chose the ForTran accent "foo[bar]" instead of the C accent "*garply" |
69 | because I felt the ForTran way was more intuitive. The C way would | |
70 | probably yield better code on most C compilers. Dean Elsner. | |
f851444e | 71 | (C style also gives deeper insight [to me] ... oh well ...) */ |
252b5132 | 72 | \f |
f851444e | 73 | void |
dd625418 KH |
74 | flonum_multip (const FLONUM_TYPE *a, const FLONUM_TYPE *b, |
75 | FLONUM_TYPE *product) | |
252b5132 | 76 | { |
f851444e NC |
77 | int size_of_a; /* 0 origin */ |
78 | int size_of_b; /* 0 origin */ | |
79 | int size_of_product; /* 0 origin */ | |
80 | int size_of_sum; /* 0 origin */ | |
81 | int extra_product_positions; /* 1 origin */ | |
252b5132 RH |
82 | unsigned long work; |
83 | unsigned long carry; | |
84 | long exponent; | |
85 | LITTLENUM_TYPE *q; | |
86 | long significant; /* TRUE when we emit a non-0 littlenum */ | |
f851444e NC |
87 | /* ForTran accent follows. */ |
88 | int P; /* Scan product low-order -> high. */ | |
252b5132 | 89 | int N; /* As in sum above. */ |
f851444e NC |
90 | int A; /* Which [] of a? */ |
91 | int B; /* Which [] of b? */ | |
252b5132 | 92 | |
f851444e NC |
93 | if ((a->sign != '-' && a->sign != '+') |
94 | || (b->sign != '-' && b->sign != '+')) | |
252b5132 | 95 | { |
f851444e | 96 | /* Got to fail somehow. Any suggestions? */ |
252b5132 RH |
97 | product->sign = 0; |
98 | return; | |
99 | } | |
100 | product->sign = (a->sign == b->sign) ? '+' : '-'; | |
101 | size_of_a = a->leader - a->low; | |
102 | size_of_b = b->leader - b->low; | |
103 | exponent = a->exponent + b->exponent; | |
104 | size_of_product = product->high - product->low; | |
105 | size_of_sum = size_of_a + size_of_b; | |
106 | extra_product_positions = size_of_product - size_of_sum; | |
107 | if (extra_product_positions < 0) | |
108 | { | |
f851444e NC |
109 | P = extra_product_positions; /* P < 0 */ |
110 | exponent -= extra_product_positions; /* Increases exponent. */ | |
252b5132 RH |
111 | } |
112 | else | |
113 | { | |
114 | P = 0; | |
115 | } | |
116 | carry = 0; | |
117 | significant = 0; | |
118 | for (N = 0; N <= size_of_sum; N++) | |
119 | { | |
120 | work = carry; | |
121 | carry = 0; | |
122 | for (A = 0; A <= N; A++) | |
123 | { | |
124 | B = N - A; | |
125 | if (A <= size_of_a && B <= size_of_b && B >= 0) | |
126 | { | |
127 | #ifdef TRACE | |
f851444e NC |
128 | printf ("a:low[%d.]=%04x b:low[%d.]=%04x work_before=%08x\n", |
129 | A, a->low[A], B, b->low[B], work); | |
252b5132 RH |
130 | #endif |
131 | /* Watch out for sign extension! Without the casts, on | |
132 | the DEC Alpha, the multiplication result is *signed* | |
133 | int, which gets sign-extended to convert to the | |
134 | unsigned long! */ | |
135 | work += (unsigned long) a->low[A] * (unsigned long) b->low[B]; | |
136 | carry += work >> LITTLENUM_NUMBER_OF_BITS; | |
137 | work &= LITTLENUM_MASK; | |
138 | #ifdef TRACE | |
139 | printf ("work=%08x carry=%04x\n", work, carry); | |
140 | #endif | |
141 | } | |
142 | } | |
143 | significant |= work; | |
144 | if (significant || P < 0) | |
145 | { | |
146 | if (P >= 0) | |
147 | { | |
148 | product->low[P] = work; | |
149 | #ifdef TRACE | |
150 | printf ("P=%d. work[p]:=%04x\n", P, work); | |
151 | #endif | |
152 | } | |
153 | P++; | |
154 | } | |
155 | else | |
156 | { | |
157 | extra_product_positions++; | |
158 | exponent++; | |
159 | } | |
160 | } | |
f851444e | 161 | /* [P]-> position # size_of_sum + 1. |
f0e652b4 | 162 | This is where 'carry' should go. */ |
252b5132 RH |
163 | #ifdef TRACE |
164 | printf ("final carry =%04x\n", carry); | |
165 | #endif | |
166 | if (carry) | |
167 | { | |
168 | if (extra_product_positions > 0) | |
f851444e | 169 | product->low[P] = carry; |
252b5132 RH |
170 | else |
171 | { | |
f851444e NC |
172 | /* No room at high order for carry littlenum. */ |
173 | /* Shift right 1 to make room for most significant littlenum. */ | |
252b5132 RH |
174 | exponent++; |
175 | P--; | |
176 | for (q = product->low + P; q >= product->low; q--) | |
177 | { | |
178 | work = *q; | |
179 | *q = carry; | |
180 | carry = work; | |
181 | } | |
182 | } | |
183 | } | |
184 | else | |
f851444e | 185 | P--; |
252b5132 RH |
186 | product->leader = product->low + P; |
187 | product->exponent = exponent; | |
188 | } |