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