1 | /*
|
---|
2 | (c) copyright 1988 by the Vrije Universiteit, Amsterdam, The Netherlands.
|
---|
3 | See the copyright notice in the ACK home directory, in the file "Copyright".
|
---|
4 | */
|
---|
5 |
|
---|
6 | /* $Header: /cvsup/minix/src/lib/ack/fphook/mul_ext.fc,v 1.1 2005/10/10 15:27:44 beng Exp $ */
|
---|
7 |
|
---|
8 | /*
|
---|
9 | ROUTINE TO MULTIPLY TWO EXTENDED FORMAT NUMBERS
|
---|
10 | */
|
---|
11 |
|
---|
12 | # include "FP_bias.h"
|
---|
13 | # include "FP_trap.h"
|
---|
14 | # include "FP_types.h"
|
---|
15 | # include "FP_shift.h"
|
---|
16 |
|
---|
17 | void
|
---|
18 | mul_ext(e1,e2)
|
---|
19 | EXTEND *e1,*e2;
|
---|
20 | {
|
---|
21 | register int i,j; /* loop control */
|
---|
22 | unsigned short mp[4]; /* multiplier */
|
---|
23 | unsigned short mc[4]; /* multipcand */
|
---|
24 | unsigned short result[8]; /* result */
|
---|
25 | register unsigned short *pres;
|
---|
26 |
|
---|
27 | /* first save the sign (XOR) */
|
---|
28 | e1->sign ^= e2->sign;
|
---|
29 |
|
---|
30 | /* compute new exponent */
|
---|
31 | e1->exp += e2->exp + 1;
|
---|
32 | /* 128 bit multiply of mantissas */
|
---|
33 |
|
---|
34 | /* assign unknown long formats */
|
---|
35 | /* to known unsigned word formats */
|
---|
36 | mp[0] = e1->m1 >> 16;
|
---|
37 | mp[1] = (unsigned short) e1->m1;
|
---|
38 | mp[2] = e1->m2 >> 16;
|
---|
39 | mp[3] = (unsigned short) e1->m2;
|
---|
40 | mc[0] = e2->m1 >> 16;
|
---|
41 | mc[1] = (unsigned short) e2->m1;
|
---|
42 | mc[2] = e2->m2 >> 16;
|
---|
43 | mc[3] = (unsigned short) e2->m2;
|
---|
44 | for (i = 8; i--;) {
|
---|
45 | result[i] = 0;
|
---|
46 | }
|
---|
47 | /*
|
---|
48 | * fill registers with their components
|
---|
49 | */
|
---|
50 | for(i=4, pres = &result[4];i--;pres--) if (mp[i]) {
|
---|
51 | unsigned short k = 0;
|
---|
52 | unsigned long mpi = mp[i];
|
---|
53 | for(j=4;j--;) {
|
---|
54 | unsigned long tmp = (unsigned long)pres[j] + k;
|
---|
55 | if (mc[j]) tmp += mpi * mc[j];
|
---|
56 | pres[j] = tmp;
|
---|
57 | k = tmp >> 16;
|
---|
58 | }
|
---|
59 | pres[-1] = k;
|
---|
60 | }
|
---|
61 | if (! (result[0] & 0x8000)) {
|
---|
62 | e1->exp--;
|
---|
63 | for (i = 0; i <= 3; i++) {
|
---|
64 | result[i] <<= 1;
|
---|
65 | if (result[i+1]&0x8000) result[i] |= 1;
|
---|
66 | }
|
---|
67 | result[4] <<= 1;
|
---|
68 | }
|
---|
69 |
|
---|
70 | /*
|
---|
71 | * combine the registers to a total
|
---|
72 | */
|
---|
73 | e1->m1 = ((unsigned long)(result[0]) << 16) + result[1];
|
---|
74 | e1->m2 = ((unsigned long)(result[2]) << 16) + result[3];
|
---|
75 | if (result[4] & 0x8000) {
|
---|
76 | if (++e1->m2 == 0)
|
---|
77 | if (++e1->m1 == 0) {
|
---|
78 | e1->m1 = NORMBIT;
|
---|
79 | e1->exp++;
|
---|
80 | }
|
---|
81 | }
|
---|
82 |
|
---|
83 | /* check for overflow */
|
---|
84 | if (e1->exp >= EXT_MAX) {
|
---|
85 | trap(EFOVFL);
|
---|
86 | /* if caught */
|
---|
87 | /* return signed infinity */
|
---|
88 | e1->exp = EXT_MAX;
|
---|
89 | infinity: e1->m1 = e1->m2 =0L;
|
---|
90 | return;
|
---|
91 | }
|
---|
92 | /* check for underflow */
|
---|
93 | if (e1->exp < EXT_MIN) {
|
---|
94 | trap(EFUNFL);
|
---|
95 | e1->exp = EXT_MIN;
|
---|
96 | goto infinity;
|
---|
97 | }
|
---|
98 | }
|
---|