annotate libdecnumber/decBasic.c @ 145:1830386684a0

gcc-9.2.0
author anatofuz
date Thu, 13 Feb 2020 11:34:05 +0900
parents 84e7813d76e9
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1 /* Common base code for the decNumber C Library.
145
1830386684a0 gcc-9.2.0
anatofuz
parents: 131
diff changeset
2 Copyright (C) 2007-2020 Free Software Foundation, Inc.
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
3 Contributed by IBM Corporation. Author Mike Cowlishaw.
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
4
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
5 This file is part of GCC.
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
6
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
7 GCC is free software; you can redistribute it and/or modify it under
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
8 the terms of the GNU General Public License as published by the Free
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
9 Software Foundation; either version 3, or (at your option) any later
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
10 version.
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
11
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
12 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
13 WARRANTY; without even the implied warranty of MERCHANTABILITY or
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
14 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
15 for more details.
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
16
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
17 Under Section 7 of GPL version 3, you are granted additional
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
18 permissions described in the GCC Runtime Library Exception, version
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
19 3.1, as published by the Free Software Foundation.
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
20
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
21 You should have received a copy of the GNU General Public License and
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
22 a copy of the GCC Runtime Library Exception along with this program;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
23 see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
24 <http://www.gnu.org/licenses/>. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
25
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
26 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
27 /* decBasic.c -- common base code for Basic decimal types */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
28 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
29 /* This module comprises code that is shared between decDouble and */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
30 /* decQuad (but not decSingle). The main arithmetic operations are */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
31 /* here (Add, Subtract, Multiply, FMA, and Division operators). */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
32 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
33 /* Unlike decNumber, parameterization takes place at compile time */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
34 /* rather than at runtime. The parameters are set in the decDouble.c */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
35 /* (etc.) files, which then include this one to produce the compiled */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
36 /* code. The functions here, therefore, are code shared between */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
37 /* multiple formats. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
38 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
39 /* This must be included after decCommon.c. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
40 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
41 /* Names here refer to decFloat rather than to decDouble, etc., and */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
42 /* the functions are in strict alphabetical order. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
43
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
44 /* The compile-time flags SINGLE, DOUBLE, and QUAD are set up in */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
45 /* decCommon.c */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
46 #if !defined(QUAD)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
47 #error decBasic.c must be included after decCommon.c
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
48 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
49 #if SINGLE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
50 #error Routines in decBasic.c are for decDouble and decQuad only
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
51 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
52
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
53 /* Private constants */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
54 #define DIVIDE 0x80000000 /* Divide operations [as flags] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
55 #define REMAINDER 0x40000000 /* .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
56 #define DIVIDEINT 0x20000000 /* .. */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
57 #define REMNEAR 0x10000000 /* .. */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
58
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
59 /* Private functions (local, used only by routines in this module) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
60 static decFloat *decDivide(decFloat *, const decFloat *,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
61 const decFloat *, decContext *, uInt);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
62 static decFloat *decCanonical(decFloat *, const decFloat *);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
63 static void decFiniteMultiply(bcdnum *, uByte *, const decFloat *,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
64 const decFloat *);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
65 static decFloat *decInfinity(decFloat *, const decFloat *);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
66 static decFloat *decInvalid(decFloat *, decContext *);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
67 static decFloat *decNaNs(decFloat *, const decFloat *, const decFloat *,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
68 decContext *);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
69 static Int decNumCompare(const decFloat *, const decFloat *, Flag);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
70 static decFloat *decToIntegral(decFloat *, const decFloat *, decContext *,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
71 enum rounding, Flag);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
72 static uInt decToInt32(const decFloat *, decContext *, enum rounding,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
73 Flag, Flag);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
74
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
75 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
76 /* decCanonical -- copy a decFloat, making canonical */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
77 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
78 /* result gets the canonicalized df */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
79 /* df is the decFloat to copy and make canonical */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
80 /* returns result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
81 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
82 /* This is exposed via decFloatCanonical for Double and Quad only. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
83 /* This works on specials, too; no error or exception is possible. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
84 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
85 static decFloat * decCanonical(decFloat *result, const decFloat *df) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
86 uInt encode, precode, dpd; /* work */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
87 uInt inword, uoff, canon; /* .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
88 Int n; /* counter (down) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
89 if (df!=result) *result=*df; /* effect copy if needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
90 if (DFISSPECIAL(result)) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
91 if (DFISINF(result)) return decInfinity(result, df); /* clean Infinity */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
92 /* is a NaN */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
93 DFWORD(result, 0)&=~ECONNANMASK; /* clear ECON except selector */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
94 if (DFISCCZERO(df)) return result; /* coefficient continuation is 0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
95 /* drop through to check payload */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
96 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
97 /* return quickly if the coefficient continuation is canonical */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
98 { /* declare block */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
99 #if DOUBLE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
100 uInt sourhi=DFWORD(df, 0);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
101 uInt sourlo=DFWORD(df, 1);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
102 if (CANONDPDOFF(sourhi, 8)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
103 && CANONDPDTWO(sourhi, sourlo, 30)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
104 && CANONDPDOFF(sourlo, 20)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
105 && CANONDPDOFF(sourlo, 10)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
106 && CANONDPDOFF(sourlo, 0)) return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
107 #elif QUAD
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
108 uInt sourhi=DFWORD(df, 0);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
109 uInt sourmh=DFWORD(df, 1);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
110 uInt sourml=DFWORD(df, 2);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
111 uInt sourlo=DFWORD(df, 3);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
112 if (CANONDPDOFF(sourhi, 4)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
113 && CANONDPDTWO(sourhi, sourmh, 26)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
114 && CANONDPDOFF(sourmh, 16)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
115 && CANONDPDOFF(sourmh, 6)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
116 && CANONDPDTWO(sourmh, sourml, 28)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
117 && CANONDPDOFF(sourml, 18)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
118 && CANONDPDOFF(sourml, 8)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
119 && CANONDPDTWO(sourml, sourlo, 30)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
120 && CANONDPDOFF(sourlo, 20)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
121 && CANONDPDOFF(sourlo, 10)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
122 && CANONDPDOFF(sourlo, 0)) return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
123 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
124 } /* block */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
125
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
126 /* Loop to repair a non-canonical coefficent, as needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
127 inword=DECWORDS-1; /* current input word */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
128 uoff=0; /* bit offset of declet */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
129 encode=DFWORD(result, inword);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
130 for (n=DECLETS-1; n>=0; n--) { /* count down declets of 10 bits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
131 dpd=encode>>uoff;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
132 uoff+=10;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
133 if (uoff>32) { /* crossed uInt boundary */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
134 inword--;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
135 encode=DFWORD(result, inword);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
136 uoff-=32;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
137 dpd|=encode<<(10-uoff); /* get pending bits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
138 }
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
139 dpd&=0x3ff; /* clear uninteresting bits */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
140 if (dpd<0x16e) continue; /* must be canonical */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
141 canon=BIN2DPD[DPD2BIN[dpd]]; /* determine canonical declet */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
142 if (canon==dpd) continue; /* have canonical declet */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
143 /* need to replace declet */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
144 if (uoff>=10) { /* all within current word */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
145 encode&=~(0x3ff<<(uoff-10)); /* clear the 10 bits ready for replace */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
146 encode|=canon<<(uoff-10); /* insert the canonical form */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
147 DFWORD(result, inword)=encode; /* .. and save */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
148 continue;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
149 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
150 /* straddled words */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
151 precode=DFWORD(result, inword+1); /* get previous */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
152 precode&=0xffffffff>>(10-uoff); /* clear top bits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
153 DFWORD(result, inword+1)=precode|(canon<<(32-(10-uoff)));
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
154 encode&=0xffffffff<<uoff; /* clear bottom bits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
155 encode|=canon>>(10-uoff); /* insert canonical */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
156 DFWORD(result, inword)=encode; /* .. and save */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
157 } /* n */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
158 return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
159 } /* decCanonical */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
160
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
161 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
162 /* decDivide -- divide operations */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
163 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
164 /* result gets the result of dividing dfl by dfr: */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
165 /* dfl is the first decFloat (lhs) */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
166 /* dfr is the second decFloat (rhs) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
167 /* set is the context */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
168 /* op is the operation selector */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
169 /* returns result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
170 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
171 /* op is one of DIVIDE, REMAINDER, DIVIDEINT, or REMNEAR. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
172 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
173 #define DIVCOUNT 0 /* 1 to instrument subtractions counter */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
174 #define DIVBASE ((uInt)BILLION) /* the base used for divide */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
175 #define DIVOPLEN DECPMAX9 /* operand length ('digits' base 10**9) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
176 #define DIVACCLEN (DIVOPLEN*3) /* accumulator length (ditto) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
177 static decFloat * decDivide(decFloat *result, const decFloat *dfl,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
178 const decFloat *dfr, decContext *set, uInt op) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
179 decFloat quotient; /* for remainders */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
180 bcdnum num; /* for final conversion */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
181 uInt acc[DIVACCLEN]; /* coefficent in base-billion .. */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
182 uInt div[DIVOPLEN]; /* divisor in base-billion .. */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
183 uInt quo[DIVOPLEN+1]; /* quotient in base-billion .. */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
184 uByte bcdacc[(DIVOPLEN+1)*9+2]; /* for quotient in BCD, +1, +1 */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
185 uInt *msua, *msud, *msuq; /* -> msu of acc, div, and quo */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
186 Int divunits, accunits; /* lengths */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
187 Int quodigits; /* digits in quotient */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
188 uInt *lsua, *lsuq; /* -> current acc and quo lsus */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
189 Int length, multiplier; /* work */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
190 uInt carry, sign; /* .. */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
191 uInt *ua, *ud, *uq; /* .. */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
192 uByte *ub; /* .. */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
193 uInt uiwork; /* for macros */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
194 uInt divtop; /* top unit of div adjusted for estimating */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
195 #if DIVCOUNT
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
196 static uInt maxcount=0; /* worst-seen subtractions count */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
197 uInt divcount=0; /* subtractions count [this divide] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
198 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
199
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
200 /* calculate sign */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
201 num.sign=(DFWORD(dfl, 0)^DFWORD(dfr, 0)) & DECFLOAT_Sign;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
202
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
203 if (DFISSPECIAL(dfl) || DFISSPECIAL(dfr)) { /* either is special? */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
204 /* NaNs are handled as usual */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
205 if (DFISNAN(dfl) || DFISNAN(dfr)) return decNaNs(result, dfl, dfr, set);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
206 /* one or two infinities */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
207 if (DFISINF(dfl)) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
208 if (DFISINF(dfr)) return decInvalid(result, set); /* Two infinities bad */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
209 if (op&(REMAINDER|REMNEAR)) return decInvalid(result, set); /* as is rem */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
210 /* Infinity/x is infinite and quiet, even if x=0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
211 DFWORD(result, 0)=num.sign;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
212 return decInfinity(result, result);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
213 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
214 /* must be x/Infinity -- remainders are lhs */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
215 if (op&(REMAINDER|REMNEAR)) return decCanonical(result, dfl);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
216 /* divides: return zero with correct sign and exponent depending */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
217 /* on op (Etiny for divide, 0 for divideInt) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
218 decFloatZero(result);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
219 if (op==DIVIDEINT) DFWORD(result, 0)|=num.sign; /* add sign */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
220 else DFWORD(result, 0)=num.sign; /* zeros the exponent, too */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
221 return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
222 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
223 /* next, handle zero operands (x/0 and 0/x) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
224 if (DFISZERO(dfr)) { /* x/0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
225 if (DFISZERO(dfl)) { /* 0/0 is undefined */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
226 decFloatZero(result);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
227 DFWORD(result, 0)=DECFLOAT_qNaN;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
228 set->status|=DEC_Division_undefined;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
229 return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
230 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
231 if (op&(REMAINDER|REMNEAR)) return decInvalid(result, set); /* bad rem */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
232 set->status|=DEC_Division_by_zero;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
233 DFWORD(result, 0)=num.sign;
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
234 return decInfinity(result, result); /* x/0 -> signed Infinity */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
235 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
236 num.exponent=GETEXPUN(dfl)-GETEXPUN(dfr); /* ideal exponent */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
237 if (DFISZERO(dfl)) { /* 0/x (x!=0) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
238 /* if divide, result is 0 with ideal exponent; divideInt has */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
239 /* exponent=0, remainders give zero with lower exponent */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
240 if (op&DIVIDEINT) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
241 decFloatZero(result);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
242 DFWORD(result, 0)|=num.sign; /* add sign */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
243 return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
244 }
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
245 if (!(op&DIVIDE)) { /* a remainder */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
246 /* exponent is the minimum of the operands */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
247 num.exponent=MINI(GETEXPUN(dfl), GETEXPUN(dfr));
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
248 /* if the result is zero the sign shall be sign of dfl */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
249 num.sign=DFWORD(dfl, 0)&DECFLOAT_Sign;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
250 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
251 bcdacc[0]=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
252 num.msd=bcdacc; /* -> 0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
253 num.lsd=bcdacc; /* .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
254 return decFinalize(result, &num, set); /* [divide may clamp exponent] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
255 } /* 0/x */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
256 /* [here, both operands are known to be finite and non-zero] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
257
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
258 /* extract the operand coefficents into 'units' which are */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
259 /* base-billion; the lhs is high-aligned in acc and the msu of both */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
260 /* acc and div is at the right-hand end of array (offset length-1); */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
261 /* the quotient can need one more unit than the operands as digits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
262 /* in it are not necessarily aligned neatly; further, the quotient */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
263 /* may not start accumulating until after the end of the initial */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
264 /* operand in acc if that is small (e.g., 1) so the accumulator */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
265 /* must have at least that number of units extra (at the ls end) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
266 GETCOEFFBILL(dfl, acc+DIVACCLEN-DIVOPLEN);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
267 GETCOEFFBILL(dfr, div);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
268 /* zero the low uInts of acc */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
269 acc[0]=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
270 acc[1]=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
271 acc[2]=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
272 acc[3]=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
273 #if DOUBLE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
274 #if DIVOPLEN!=2
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
275 #error Unexpected Double DIVOPLEN
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
276 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
277 #elif QUAD
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
278 acc[4]=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
279 acc[5]=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
280 acc[6]=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
281 acc[7]=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
282 #if DIVOPLEN!=4
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
283 #error Unexpected Quad DIVOPLEN
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
284 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
285 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
286
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
287 /* set msu and lsu pointers */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
288 msua=acc+DIVACCLEN-1; /* [leading zeros removed below] */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
289 msuq=quo+DIVOPLEN;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
290 /*[loop for div will terminate because operands are non-zero] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
291 for (msud=div+DIVOPLEN-1; *msud==0;) msud--;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
292 /* the initial least-significant unit of acc is set so acc appears */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
293 /* to have the same length as div. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
294 /* This moves one position towards the least possible for each */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
295 /* iteration */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
296 divunits=(Int)(msud-div+1); /* precalculate */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
297 lsua=msua-divunits+1; /* initial working lsu of acc */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
298 lsuq=msuq; /* and of quo */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
299
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
300 /* set up the estimator for the multiplier; this is the msu of div, */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
301 /* plus two bits from the unit below (if any) rounded up by one if */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
302 /* there are any non-zero bits or units below that [the extra two */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
303 /* bits makes for a much better estimate when the top unit is small] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
304 divtop=*msud<<2;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
305 if (divunits>1) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
306 uInt *um=msud-1;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
307 uInt d=*um;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
308 if (d>=750000000) {divtop+=3; d-=750000000;}
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
309 else if (d>=500000000) {divtop+=2; d-=500000000;}
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
310 else if (d>=250000000) {divtop++; d-=250000000;}
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
311 if (d) divtop++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
312 else for (um--; um>=div; um--) if (*um) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
313 divtop++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
314 break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
315 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
316 } /* >1 unit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
317
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
318 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
319 {Int i;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
320 printf("----- div=");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
321 for (i=divunits-1; i>=0; i--) printf("%09ld ", (LI)div[i]);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
322 printf("\n");}
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
323 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
324
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
325 /* now collect up to DECPMAX+1 digits in the quotient (this may */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
326 /* need OPLEN+1 uInts if unaligned) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
327 quodigits=0; /* no digits yet */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
328 for (;; lsua--) { /* outer loop -- each input position */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
329 #if DECCHECK
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
330 if (lsua<acc) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
331 printf("Acc underrun...\n");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
332 break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
333 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
334 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
335 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
336 printf("Outer: quodigits=%ld acc=", (LI)quodigits);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
337 for (ua=msua; ua>=lsua; ua--) printf("%09ld ", (LI)*ua);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
338 printf("\n");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
339 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
340 *lsuq=0; /* default unit result is 0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
341 for (;;) { /* inner loop -- calculate quotient unit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
342 /* strip leading zero units from acc (either there initially or */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
343 /* from subtraction below); this may strip all if exactly 0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
344 for (; *msua==0 && msua>=lsua;) msua--;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
345 accunits=(Int)(msua-lsua+1); /* [maybe 0] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
346 /* subtraction is only necessary and possible if there are as */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
347 /* least as many units remaining in acc for this iteration as */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
348 /* there are in div */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
349 if (accunits<divunits) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
350 if (accunits==0) msua++; /* restore */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
351 break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
352 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
353
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
354 /* If acc is longer than div then subtraction is definitely */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
355 /* possible (as msu of both is non-zero), but if they are the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
356 /* same length a comparison is needed. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
357 /* If a subtraction is needed then a good estimate of the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
358 /* multiplier for the subtraction is also needed in order to */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
359 /* minimise the iterations of this inner loop because the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
360 /* subtractions needed dominate division performance. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
361 if (accunits==divunits) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
362 /* compare the high divunits of acc and div: */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
363 /* acc<div: this quotient unit is unchanged; subtraction */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
364 /* will be possible on the next iteration */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
365 /* acc==div: quotient gains 1, set acc=0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
366 /* acc>div: subtraction necessary at this position */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
367 for (ud=msud, ua=msua; ud>div; ud--, ua--) if (*ud!=*ua) break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
368 /* [now at first mismatch or lsu] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
369 if (*ud>*ua) break; /* next time... */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
370 if (*ud==*ua) { /* all compared equal */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
371 *lsuq+=1; /* increment result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
372 msua=lsua; /* collapse acc units */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
373 *msua=0; /* .. to a zero */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
374 break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
375 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
376
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
377 /* subtraction necessary; estimate multiplier [see above] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
378 /* if both *msud and *msua are small it is cost-effective to */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
379 /* bring in part of the following units (if any) to get a */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
380 /* better estimate (assume some other non-zero in div) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
381 #define DIVLO 1000000U
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
382 #define DIVHI (DIVBASE/DIVLO)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
383 #if DECUSE64
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
384 if (divunits>1) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
385 /* there cannot be a *(msud-2) for DECDOUBLE so next is */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
386 /* an exact calculation unless DECQUAD (which needs to */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
387 /* assume bits out there if divunits>2) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
388 uLong mul=(uLong)*msua * DIVBASE + *(msua-1);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
389 uLong div=(uLong)*msud * DIVBASE + *(msud-1);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
390 #if QUAD
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
391 if (divunits>2) div++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
392 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
393 mul/=div;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
394 multiplier=(Int)mul;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
395 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
396 else multiplier=*msua/(*msud);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
397 #else
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
398 if (divunits>1 && *msua<DIVLO && *msud<DIVLO) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
399 multiplier=(*msua*DIVHI + *(msua-1)/DIVLO)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
400 /(*msud*DIVHI + *(msud-1)/DIVLO +1);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
401 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
402 else multiplier=(*msua<<2)/divtop;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
403 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
404 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
405 else { /* accunits>divunits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
406 /* msud is one unit 'lower' than msua, so estimate differently */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
407 #if DECUSE64
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
408 uLong mul;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
409 /* as before, bring in extra digits if possible */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
410 if (divunits>1 && *msua<DIVLO && *msud<DIVLO) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
411 mul=((uLong)*msua * DIVHI * DIVBASE) + *(msua-1) * DIVHI
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
412 + *(msua-2)/DIVLO;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
413 mul/=(*msud*DIVHI + *(msud-1)/DIVLO +1);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
414 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
415 else if (divunits==1) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
416 mul=(uLong)*msua * DIVBASE + *(msua-1);
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
417 mul/=*msud; /* no more to the right */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
418 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
419 else {
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
420 mul=(uLong)(*msua) * (uInt)(DIVBASE<<2)
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
421 + (*(msua-1)<<2);
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
422 mul/=divtop; /* [divtop already allows for sticky bits] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
423 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
424 multiplier=(Int)mul;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
425 #else
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
426 multiplier=*msua * ((DIVBASE<<2)/divtop);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
427 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
428 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
429 if (multiplier==0) multiplier=1; /* marginal case */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
430 *lsuq+=multiplier;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
431
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
432 #if DIVCOUNT
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
433 /* printf("Multiplier: %ld\n", (LI)multiplier); */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
434 divcount++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
435 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
436
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
437 /* Carry out the subtraction acc-(div*multiplier); for each */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
438 /* unit in div, do the multiply, split to units (see */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
439 /* decFloatMultiply for the algorithm), and subtract from acc */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
440 #define DIVMAGIC 2305843009U /* 2**61/10**9 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
441 #define DIVSHIFTA 29
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
442 #define DIVSHIFTB 32
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
443 carry=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
444 for (ud=div, ua=lsua; ud<=msud; ud++, ua++) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
445 uInt lo, hop;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
446 #if DECUSE64
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
447 uLong sub=(uLong)multiplier*(*ud)+carry;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
448 if (sub<DIVBASE) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
449 carry=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
450 lo=(uInt)sub;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
451 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
452 else {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
453 hop=(uInt)(sub>>DIVSHIFTA);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
454 carry=(uInt)(((uLong)hop*DIVMAGIC)>>DIVSHIFTB);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
455 /* the estimate is now in hi; now calculate sub-hi*10**9 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
456 /* to get the remainder (which will be <DIVBASE)) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
457 lo=(uInt)sub;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
458 lo-=carry*DIVBASE; /* low word of result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
459 if (lo>=DIVBASE) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
460 lo-=DIVBASE; /* correct by +1 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
461 carry++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
462 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
463 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
464 #else /* 32-bit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
465 uInt hi;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
466 /* calculate multiplier*(*ud) into hi and lo */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
467 LONGMUL32HI(hi, *ud, multiplier); /* get the high word */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
468 lo=multiplier*(*ud); /* .. and the low */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
469 lo+=carry; /* add the old hi */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
470 carry=hi+(lo<carry); /* .. with any carry */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
471 if (carry || lo>=DIVBASE) { /* split is needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
472 hop=(carry<<3)+(lo>>DIVSHIFTA); /* hi:lo/2**29 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
473 LONGMUL32HI(carry, hop, DIVMAGIC); /* only need the high word */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
474 /* [DIVSHIFTB is 32, so carry can be used directly] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
475 /* the estimate is now in carry; now calculate hi:lo-est*10**9; */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
476 /* happily the top word of the result is irrelevant because it */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
477 /* will always be zero so this needs only one multiplication */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
478 lo-=(carry*DIVBASE);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
479 /* the correction here will be at most +1; do it */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
480 if (lo>=DIVBASE) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
481 lo-=DIVBASE;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
482 carry++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
483 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
484 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
485 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
486 if (lo>*ua) { /* borrow needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
487 *ua+=DIVBASE;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
488 carry++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
489 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
490 *ua-=lo;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
491 } /* ud loop */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
492 if (carry) *ua-=carry; /* accdigits>divdigits [cannot borrow] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
493 } /* inner loop */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
494
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
495 /* the outer loop terminates when there is either an exact result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
496 /* or enough digits; first update the quotient digit count and */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
497 /* pointer (if any significant digits) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
498 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
499 if (*lsuq || quodigits) printf("*lsuq=%09ld\n", (LI)*lsuq);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
500 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
501 if (quodigits) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
502 quodigits+=9; /* had leading unit earlier */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
503 lsuq--;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
504 if (quodigits>DECPMAX+1) break; /* have enough */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
505 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
506 else if (*lsuq) { /* first quotient digits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
507 const uInt *pow;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
508 for (pow=DECPOWERS; *lsuq>=*pow; pow++) quodigits++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
509 lsuq--;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
510 /* [cannot have >DECPMAX+1 on first unit] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
511 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
512
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
513 if (*msua!=0) continue; /* not an exact result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
514 /* acc is zero iff used all of original units and zero down to lsua */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
515 /* (must also continue to original lsu for correct quotient length) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
516 if (lsua>acc+DIVACCLEN-DIVOPLEN) continue;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
517 for (; msua>lsua && *msua==0;) msua--;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
518 if (*msua==0 && msua==lsua) break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
519 } /* outer loop */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
520
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
521 /* all of the original operand in acc has been covered at this point */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
522 /* quotient now has at least DECPMAX+2 digits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
523 /* *msua is now non-0 if inexact and sticky bits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
524 /* lsuq is one below the last uint of the quotient */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
525 lsuq++; /* set -> true lsu of quo */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
526 if (*msua) *lsuq|=1; /* apply sticky bit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
527
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
528 /* quo now holds the (unrounded) quotient in base-billion; one */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
529 /* base-billion 'digit' per uInt. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
530 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
531 printf("DivQuo:");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
532 for (uq=msuq; uq>=lsuq; uq--) printf(" %09ld", (LI)*uq);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
533 printf("\n");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
534 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
535
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
536 /* Now convert to BCD for rounding and cleanup, starting from the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
537 /* most significant end [offset by one into bcdacc to leave room */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
538 /* for a possible carry digit if rounding for REMNEAR is needed] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
539 for (uq=msuq, ub=bcdacc+1; uq>=lsuq; uq--, ub+=9) {
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
540 uInt top, mid, rem; /* work */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
541 if (*uq==0) { /* no split needed */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
542 UBFROMUI(ub, 0); /* clear 9 BCD8s */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
543 UBFROMUI(ub+4, 0); /* .. */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
544 *(ub+8)=0; /* .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
545 continue;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
546 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
547 /* *uq is non-zero -- split the base-billion digit into */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
548 /* hi, mid, and low three-digits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
549 #define divsplit9 1000000 /* divisor */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
550 #define divsplit6 1000 /* divisor */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
551 /* The splitting is done by simple divides and remainders, */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
552 /* assuming the compiler will optimize these [GCC does] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
553 top=*uq/divsplit9;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
554 rem=*uq%divsplit9;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
555 mid=rem/divsplit6;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
556 rem=rem%divsplit6;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
557 /* lay out the nine BCD digits (plus one unwanted byte) */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
558 UBFROMUI(ub, UBTOUI(&BIN2BCD8[top*4]));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
559 UBFROMUI(ub+3, UBTOUI(&BIN2BCD8[mid*4]));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
560 UBFROMUI(ub+6, UBTOUI(&BIN2BCD8[rem*4]));
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
561 } /* BCD conversion loop */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
562 ub--; /* -> lsu */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
563
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
564 /* complete the bcdnum; quodigits is correct, so the position of */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
565 /* the first non-zero is known */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
566 num.msd=bcdacc+1+(msuq-lsuq+1)*9-quodigits;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
567 num.lsd=ub;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
568
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
569 /* make exponent adjustments, etc */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
570 if (lsua<acc+DIVACCLEN-DIVOPLEN) { /* used extra digits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
571 num.exponent-=(Int)((acc+DIVACCLEN-DIVOPLEN-lsua)*9);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
572 /* if the result was exact then there may be up to 8 extra */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
573 /* trailing zeros in the overflowed quotient final unit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
574 if (*msua==0) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
575 for (; *ub==0;) ub--; /* drop zeros */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
576 num.exponent+=(Int)(num.lsd-ub); /* and adjust exponent */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
577 num.lsd=ub;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
578 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
579 } /* adjustment needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
580
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
581 #if DIVCOUNT
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
582 if (divcount>maxcount) { /* new high-water nark */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
583 maxcount=divcount;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
584 printf("DivNewMaxCount: %ld\n", (LI)maxcount);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
585 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
586 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
587
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
588 if (op&DIVIDE) return decFinalize(result, &num, set); /* all done */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
589
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
590 /* Is DIVIDEINT or a remainder; there is more to do -- first form */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
591 /* the integer (this is done 'after the fact', unlike as in */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
592 /* decNumber, so as not to tax DIVIDE) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
593
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
594 /* The first non-zero digit will be in the first 9 digits, known */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
595 /* from quodigits and num.msd, so there is always space for DECPMAX */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
596 /* digits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
597
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
598 length=(Int)(num.lsd-num.msd+1);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
599 /*printf("Length exp: %ld %ld\n", (LI)length, (LI)num.exponent); */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
600
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
601 if (length+num.exponent>DECPMAX) { /* cannot fit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
602 decFloatZero(result);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
603 DFWORD(result, 0)=DECFLOAT_qNaN;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
604 set->status|=DEC_Division_impossible;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
605 return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
606 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
607
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
608 if (num.exponent>=0) { /* already an int, or need pad zeros */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
609 for (ub=num.lsd+1; ub<=num.lsd+num.exponent; ub++) *ub=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
610 num.lsd+=num.exponent;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
611 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
612 else { /* too long: round or truncate needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
613 Int drop=-num.exponent;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
614 if (!(op&REMNEAR)) { /* simple truncate */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
615 num.lsd-=drop;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
616 if (num.lsd<num.msd) { /* truncated all */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
617 num.lsd=num.msd; /* make 0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
618 *num.lsd=0; /* .. [sign still relevant] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
619 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
620 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
621 else { /* round to nearest even [sigh] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
622 /* round-to-nearest, in-place; msd is at or to right of bcdacc+1 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
623 /* (this is a special case of Quantize -- q.v. for commentary) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
624 uByte *roundat; /* -> re-round digit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
625 uByte reround; /* reround value */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
626 *(num.msd-1)=0; /* in case of left carry, or make 0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
627 if (drop<length) roundat=num.lsd-drop+1;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
628 else if (drop==length) roundat=num.msd;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
629 else roundat=num.msd-1; /* [-> 0] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
630 reround=*roundat;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
631 for (ub=roundat+1; ub<=num.lsd; ub++) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
632 if (*ub!=0) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
633 reround=DECSTICKYTAB[reround];
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
634 break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
635 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
636 } /* check stickies */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
637 if (roundat>num.msd) num.lsd=roundat-1;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
638 else {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
639 num.msd--; /* use the 0 .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
640 num.lsd=num.msd; /* .. at the new MSD place */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
641 }
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
642 if (reround!=0) { /* discarding non-zero */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
643 uInt bump=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
644 /* rounding is DEC_ROUND_HALF_EVEN always */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
645 if (reround>5) bump=1; /* >0.5 goes up */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
646 else if (reround==5) /* exactly 0.5000 .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
647 bump=*(num.lsd) & 0x01; /* .. up iff [new] lsd is odd */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
648 if (bump!=0) { /* need increment */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
649 /* increment the coefficient; this might end up with 1000... */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
650 ub=num.lsd;
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
651 for (; UBTOUI(ub-3)==0x09090909; ub-=4) UBFROMUI(ub-3, 0);
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
652 for (; *ub==9; ub--) *ub=0; /* at most 3 more */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
653 *ub+=1;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
654 if (ub<num.msd) num.msd--; /* carried */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
655 } /* bump needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
656 } /* reround!=0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
657 } /* remnear */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
658 } /* round or truncate needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
659 num.exponent=0; /* all paths */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
660 /*decShowNum(&num, "int"); */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
661
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
662 if (op&DIVIDEINT) return decFinalize(result, &num, set); /* all done */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
663
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
664 /* Have a remainder to calculate */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
665 decFinalize(&quotient, &num, set); /* lay out the integer so far */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
666 DFWORD(&quotient, 0)^=DECFLOAT_Sign; /* negate it */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
667 sign=DFWORD(dfl, 0); /* save sign of dfl */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
668 decFloatFMA(result, &quotient, dfr, dfl, set);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
669 if (!DFISZERO(result)) return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
670 /* if the result is zero the sign shall be sign of dfl */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
671 DFWORD(&quotient, 0)=sign; /* construct decFloat of sign */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
672 return decFloatCopySign(result, result, &quotient);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
673 } /* decDivide */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
674
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
675 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
676 /* decFiniteMultiply -- multiply two finite decFloats */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
677 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
678 /* num gets the result of multiplying dfl and dfr */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
679 /* bcdacc .. with the coefficient in this array */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
680 /* dfl is the first decFloat (lhs) */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
681 /* dfr is the second decFloat (rhs) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
682 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
683 /* This effects the multiplication of two decFloats, both known to be */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
684 /* finite, leaving the result in a bcdnum ready for decFinalize (for */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
685 /* use in Multiply) or in a following addition (FMA). */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
686 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
687 /* bcdacc must have space for at least DECPMAX9*18+1 bytes. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
688 /* No error is possible and no status is set. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
689 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
690 /* This routine has two separate implementations of the core */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
691 /* multiplication; both using base-billion. One uses only 32-bit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
692 /* variables (Ints and uInts) or smaller; the other uses uLongs (for */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
693 /* multiplication and addition only). Both implementations cover */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
694 /* both arithmetic sizes (DOUBLE and QUAD) in order to allow timing */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
695 /* comparisons. In any one compilation only one implementation for */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
696 /* each size can be used, and if DECUSE64 is 0 then use of the 32-bit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
697 /* version is forced. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
698 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
699 /* Historical note: an earlier version of this code also supported the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
700 /* 256-bit format and has been preserved. That is somewhat trickier */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
701 /* during lazy carry splitting because the initial quotient estimate */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
702 /* (est) can exceed 32 bits. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
703
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
704 #define MULTBASE ((uInt)BILLION) /* the base used for multiply */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
705 #define MULOPLEN DECPMAX9 /* operand length ('digits' base 10**9) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
706 #define MULACCLEN (MULOPLEN*2) /* accumulator length (ditto) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
707 #define LEADZEROS (MULACCLEN*9 - DECPMAX*2) /* leading zeros always */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
708
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
709 /* Assertions: exponent not too large and MULACCLEN is a multiple of 4 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
710 #if DECEMAXD>9
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
711 #error Exponent may overflow when doubled for Multiply
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
712 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
713 #if MULACCLEN!=(MULACCLEN/4)*4
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
714 /* This assumption is used below only for initialization */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
715 #error MULACCLEN is not a multiple of 4
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
716 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
717
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
718 static void decFiniteMultiply(bcdnum *num, uByte *bcdacc,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
719 const decFloat *dfl, const decFloat *dfr) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
720 uInt bufl[MULOPLEN]; /* left coefficient (base-billion) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
721 uInt bufr[MULOPLEN]; /* right coefficient (base-billion) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
722 uInt *ui, *uj; /* work */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
723 uByte *ub; /* .. */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
724 uInt uiwork; /* for macros */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
725
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
726 #if DECUSE64
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
727 uLong accl[MULACCLEN]; /* lazy accumulator (base-billion+) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
728 uLong *pl; /* work -> lazy accumulator */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
729 uInt acc[MULACCLEN]; /* coefficent in base-billion .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
730 #else
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
731 uInt acc[MULACCLEN*2]; /* accumulator in base-billion .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
732 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
733 uInt *pa; /* work -> accumulator */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
734 /*printf("Base10**9: OpLen=%d MulAcclen=%d\n", OPLEN, MULACCLEN); */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
735
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
736 /* Calculate sign and exponent */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
737 num->sign=(DFWORD(dfl, 0)^DFWORD(dfr, 0)) & DECFLOAT_Sign;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
738 num->exponent=GETEXPUN(dfl)+GETEXPUN(dfr); /* [see assertion above] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
739
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
740 /* Extract the coefficients and prepare the accumulator */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
741 /* the coefficients of the operands are decoded into base-billion */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
742 /* numbers in uInt arrays (bufl and bufr, LSD at offset 0) of the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
743 /* appropriate size. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
744 GETCOEFFBILL(dfl, bufl);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
745 GETCOEFFBILL(dfr, bufr);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
746 #if DECTRACE && 0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
747 printf("CoeffbL:");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
748 for (ui=bufl+MULOPLEN-1; ui>=bufl; ui--) printf(" %08lx", (LI)*ui);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
749 printf("\n");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
750 printf("CoeffbR:");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
751 for (uj=bufr+MULOPLEN-1; uj>=bufr; uj--) printf(" %08lx", (LI)*uj);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
752 printf("\n");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
753 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
754
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
755 /* start the 64-bit/32-bit differing paths... */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
756 #if DECUSE64
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
757
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
758 /* zero the accumulator */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
759 #if MULACCLEN==4
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
760 accl[0]=0; accl[1]=0; accl[2]=0; accl[3]=0;
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
761 #else /* use a loop */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
762 /* MULACCLEN is a multiple of four, asserted above */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
763 for (pl=accl; pl<accl+MULACCLEN; pl+=4) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
764 *pl=0; *(pl+1)=0; *(pl+2)=0; *(pl+3)=0;/* [reduce overhead] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
765 } /* pl */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
766 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
767
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
768 /* Effect the multiplication */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
769 /* The multiplcation proceeds using MFC's lazy-carry resolution */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
770 /* algorithm from decNumber. First, the multiplication is */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
771 /* effected, allowing accumulation of the partial products (which */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
772 /* are in base-billion at each column position) into 64 bits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
773 /* without resolving back to base=billion after each addition. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
774 /* These 64-bit numbers (which may contain up to 19 decimal digits) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
775 /* are then split using the Clark & Cowlishaw algorithm (see below). */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
776 /* [Testing for 0 in the inner loop is not really a 'win'] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
777 for (ui=bufr; ui<bufr+MULOPLEN; ui++) { /* over each item in rhs */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
778 if (*ui==0) continue; /* product cannot affect result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
779 pl=accl+(ui-bufr); /* where to add the lhs */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
780 for (uj=bufl; uj<bufl+MULOPLEN; uj++, pl++) { /* over each item in lhs */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
781 /* if (*uj==0) continue; // product cannot affect result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
782 *pl+=((uLong)*ui)*(*uj);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
783 } /* uj */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
784 } /* ui */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
785
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
786 /* The 64-bit carries must now be resolved; this means that a */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
787 /* quotient/remainder has to be calculated for base-billion (1E+9). */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
788 /* For this, Clark & Cowlishaw's quotient estimation approach (also */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
789 /* used in decNumber) is needed, because 64-bit divide is generally */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
790 /* extremely slow on 32-bit machines, and may be slower than this */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
791 /* approach even on 64-bit machines. This algorithm splits X */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
792 /* using: */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
793 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
794 /* magic=2**(A+B)/1E+9; // 'magic number' */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
795 /* hop=X/2**A; // high order part of X (by shift) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
796 /* est=magic*hop/2**B // quotient estimate (may be low by 1) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
797 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
798 /* A and B are quite constrained; hop and magic must fit in 32 bits, */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
799 /* and 2**(A+B) must be as large as possible (which is 2**61 if */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
800 /* magic is to fit). Further, maxX increases with the length of */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
801 /* the operands (and hence the number of partial products */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
802 /* accumulated); maxX is OPLEN*(10**18), which is up to 19 digits. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
803 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
804 /* It can be shown that when OPLEN is 2 then the maximum error in */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
805 /* the estimated quotient is <1, but for larger maximum x the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
806 /* maximum error is above 1 so a correction that is >1 may be */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
807 /* needed. Values of A and B are chosen to satisfy the constraints */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
808 /* just mentioned while minimizing the maximum error (and hence the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
809 /* maximum correction), as shown in the following table: */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
810 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
811 /* Type OPLEN A B maxX maxError maxCorrection */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
812 /* --------------------------------------------------------- */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
813 /* DOUBLE 2 29 32 <2*10**18 0.63 1 */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
814 /* QUAD 4 30 31 <4*10**18 1.17 2 */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
815 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
816 /* In the OPLEN==2 case there is most choice, but the value for B */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
817 /* of 32 has a big advantage as then the calculation of the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
818 /* estimate requires no shifting; the compiler can extract the high */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
819 /* word directly after multiplying magic*hop. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
820 #define MULMAGIC 2305843009U /* 2**61/10**9 [both cases] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
821 #if DOUBLE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
822 #define MULSHIFTA 29
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
823 #define MULSHIFTB 32
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
824 #elif QUAD
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
825 #define MULSHIFTA 30
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
826 #define MULSHIFTB 31
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
827 #else
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
828 #error Unexpected type
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
829 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
830
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
831 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
832 printf("MulAccl:");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
833 for (pl=accl+MULACCLEN-1; pl>=accl; pl--)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
834 printf(" %08lx:%08lx", (LI)(*pl>>32), (LI)(*pl&0xffffffff));
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
835 printf("\n");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
836 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
837
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
838 for (pl=accl, pa=acc; pl<accl+MULACCLEN; pl++, pa++) { /* each column position */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
839 uInt lo, hop; /* work */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
840 uInt est; /* cannot exceed 4E+9 */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
841 if (*pl>=MULTBASE) {
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
842 /* *pl holds a binary number which needs to be split */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
843 hop=(uInt)(*pl>>MULSHIFTA);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
844 est=(uInt)(((uLong)hop*MULMAGIC)>>MULSHIFTB);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
845 /* the estimate is now in est; now calculate hi:lo-est*10**9; */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
846 /* happily the top word of the result is irrelevant because it */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
847 /* will always be zero so this needs only one multiplication */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
848 lo=(uInt)(*pl-((uLong)est*MULTBASE)); /* low word of result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
849 /* If QUAD, the correction here could be +2 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
850 if (lo>=MULTBASE) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
851 lo-=MULTBASE; /* correct by +1 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
852 est++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
853 #if QUAD
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
854 /* may need to correct by +2 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
855 if (lo>=MULTBASE) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
856 lo-=MULTBASE;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
857 est++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
858 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
859 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
860 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
861 /* finally place lo as the new coefficient 'digit' and add est to */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
862 /* the next place up [this is safe because this path is never */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
863 /* taken on the final iteration as *pl will fit] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
864 *pa=lo;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
865 *(pl+1)+=est;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
866 } /* *pl needed split */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
867 else { /* *pl<MULTBASE */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
868 *pa=(uInt)*pl; /* just copy across */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
869 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
870 } /* pl loop */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
871
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
872 #else /* 32-bit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
873 for (pa=acc;; pa+=4) { /* zero the accumulator */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
874 *pa=0; *(pa+1)=0; *(pa+2)=0; *(pa+3)=0; /* [reduce overhead] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
875 if (pa==acc+MULACCLEN*2-4) break; /* multiple of 4 asserted */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
876 } /* pa */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
877
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
878 /* Effect the multiplication */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
879 /* uLongs are not available (and in particular, there is no uLong */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
880 /* divide) but it is still possible to use MFC's lazy-carry */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
881 /* resolution algorithm from decNumber. First, the multiplication */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
882 /* is effected, allowing accumulation of the partial products */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
883 /* (which are in base-billion at each column position) into 64 bits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
884 /* [with the high-order 32 bits in each position being held at */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
885 /* offset +ACCLEN from the low-order 32 bits in the accumulator]. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
886 /* These 64-bit numbers (which may contain up to 19 decimal digits) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
887 /* are then split using the Clark & Cowlishaw algorithm (see */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
888 /* below). */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
889 for (ui=bufr;; ui++) { /* over each item in rhs */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
890 uInt hi, lo; /* words of exact multiply result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
891 pa=acc+(ui-bufr); /* where to add the lhs */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
892 for (uj=bufl;; uj++, pa++) { /* over each item in lhs */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
893 LONGMUL32HI(hi, *ui, *uj); /* calculate product of digits */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
894 lo=(*ui)*(*uj); /* .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
895 *pa+=lo; /* accumulate low bits and .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
896 *(pa+MULACCLEN)+=hi+(*pa<lo); /* .. high bits with any carry */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
897 if (uj==bufl+MULOPLEN-1) break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
898 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
899 if (ui==bufr+MULOPLEN-1) break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
900 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
901
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
902 /* The 64-bit carries must now be resolved; this means that a */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
903 /* quotient/remainder has to be calculated for base-billion (1E+9). */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
904 /* For this, Clark & Cowlishaw's quotient estimation approach (also */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
905 /* used in decNumber) is needed, because 64-bit divide is generally */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
906 /* extremely slow on 32-bit machines. This algorithm splits X */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
907 /* using: */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
908 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
909 /* magic=2**(A+B)/1E+9; // 'magic number' */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
910 /* hop=X/2**A; // high order part of X (by shift) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
911 /* est=magic*hop/2**B // quotient estimate (may be low by 1) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
912 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
913 /* A and B are quite constrained; hop and magic must fit in 32 bits, */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
914 /* and 2**(A+B) must be as large as possible (which is 2**61 if */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
915 /* magic is to fit). Further, maxX increases with the length of */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
916 /* the operands (and hence the number of partial products */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
917 /* accumulated); maxX is OPLEN*(10**18), which is up to 19 digits. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
918 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
919 /* It can be shown that when OPLEN is 2 then the maximum error in */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
920 /* the estimated quotient is <1, but for larger maximum x the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
921 /* maximum error is above 1 so a correction that is >1 may be */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
922 /* needed. Values of A and B are chosen to satisfy the constraints */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
923 /* just mentioned while minimizing the maximum error (and hence the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
924 /* maximum correction), as shown in the following table: */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
925 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
926 /* Type OPLEN A B maxX maxError maxCorrection */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
927 /* --------------------------------------------------------- */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
928 /* DOUBLE 2 29 32 <2*10**18 0.63 1 */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
929 /* QUAD 4 30 31 <4*10**18 1.17 2 */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
930 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
931 /* In the OPLEN==2 case there is most choice, but the value for B */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
932 /* of 32 has a big advantage as then the calculation of the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
933 /* estimate requires no shifting; the high word is simply */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
934 /* calculated from multiplying magic*hop. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
935 #define MULMAGIC 2305843009U /* 2**61/10**9 [both cases] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
936 #if DOUBLE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
937 #define MULSHIFTA 29
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
938 #define MULSHIFTB 32
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
939 #elif QUAD
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
940 #define MULSHIFTA 30
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
941 #define MULSHIFTB 31
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
942 #else
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
943 #error Unexpected type
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
944 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
945
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
946 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
947 printf("MulHiLo:");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
948 for (pa=acc+MULACCLEN-1; pa>=acc; pa--)
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
949 printf(" %08lx:%08lx", (LI)*(pa+MULACCLEN), (LI)*pa);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
950 printf("\n");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
951 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
952
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
953 for (pa=acc;; pa++) { /* each low uInt */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
954 uInt hi, lo; /* words of exact multiply result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
955 uInt hop, estlo; /* work */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
956 #if QUAD
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
957 uInt esthi; /* .. */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
958 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
959
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
960 lo=*pa;
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
961 hi=*(pa+MULACCLEN); /* top 32 bits */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
962 /* hi and lo now hold a binary number which needs to be split */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
963
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
964 #if DOUBLE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
965 hop=(hi<<3)+(lo>>MULSHIFTA); /* hi:lo/2**29 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
966 LONGMUL32HI(estlo, hop, MULMAGIC);/* only need the high word */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
967 /* [MULSHIFTB is 32, so estlo can be used directly] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
968 /* the estimate is now in estlo; now calculate hi:lo-est*10**9; */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
969 /* happily the top word of the result is irrelevant because it */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
970 /* will always be zero so this needs only one multiplication */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
971 lo-=(estlo*MULTBASE);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
972 /* esthi=0; // high word is ignored below */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
973 /* the correction here will be at most +1; do it */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
974 if (lo>=MULTBASE) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
975 lo-=MULTBASE;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
976 estlo++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
977 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
978 #elif QUAD
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
979 hop=(hi<<2)+(lo>>MULSHIFTA); /* hi:lo/2**30 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
980 LONGMUL32HI(esthi, hop, MULMAGIC);/* shift will be 31 .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
981 estlo=hop*MULMAGIC; /* .. so low word needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
982 estlo=(esthi<<1)+(estlo>>MULSHIFTB); /* [just the top bit] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
983 /* esthi=0; // high word is ignored below */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
984 lo-=(estlo*MULTBASE); /* as above */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
985 /* the correction here could be +1 or +2 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
986 if (lo>=MULTBASE) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
987 lo-=MULTBASE;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
988 estlo++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
989 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
990 if (lo>=MULTBASE) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
991 lo-=MULTBASE;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
992 estlo++;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
993 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
994 #else
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
995 #error Unexpected type
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
996 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
997
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
998 /* finally place lo as the new accumulator digit and add est to */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
999 /* the next place up; this latter add could cause a carry of 1 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1000 /* to the high word of the next place */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1001 *pa=lo;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1002 *(pa+1)+=estlo;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1003 /* esthi is always 0 for DOUBLE and QUAD so this is skipped */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1004 /* *(pa+1+MULACCLEN)+=esthi; */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1005 if (*(pa+1)<estlo) *(pa+1+MULACCLEN)+=1; /* carry */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1006 if (pa==acc+MULACCLEN-2) break; /* [MULACCLEN-1 will never need split] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1007 } /* pa loop */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1008 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1009
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1010 /* At this point, whether using the 64-bit or the 32-bit paths, the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1011 /* accumulator now holds the (unrounded) result in base-billion; */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1012 /* one base-billion 'digit' per uInt. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1013 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1014 printf("MultAcc:");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1015 for (pa=acc+MULACCLEN-1; pa>=acc; pa--) printf(" %09ld", (LI)*pa);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1016 printf("\n");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1017 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1018
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1019 /* Now convert to BCD for rounding and cleanup, starting from the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1020 /* most significant end */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1021 pa=acc+MULACCLEN-1;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1022 if (*pa!=0) num->msd=bcdacc+LEADZEROS;/* drop known lead zeros */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1023 else { /* >=1 word of leading zeros */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1024 num->msd=bcdacc; /* known leading zeros are gone */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1025 pa--; /* skip first word .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1026 for (; *pa==0; pa--) if (pa==acc) break; /* .. and any more leading 0s */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1027 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1028 for (ub=bcdacc;; pa--, ub+=9) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1029 if (*pa!=0) { /* split(s) needed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1030 uInt top, mid, rem; /* work */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1031 /* *pa is non-zero -- split the base-billion acc digit into */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1032 /* hi, mid, and low three-digits */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1033 #define mulsplit9 1000000 /* divisor */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1034 #define mulsplit6 1000 /* divisor */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1035 /* The splitting is done by simple divides and remainders, */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1036 /* assuming the compiler will optimize these where useful */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1037 /* [GCC does] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1038 top=*pa/mulsplit9;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1039 rem=*pa%mulsplit9;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1040 mid=rem/mulsplit6;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1041 rem=rem%mulsplit6;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1042 /* lay out the nine BCD digits (plus one unwanted byte) */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1043 UBFROMUI(ub, UBTOUI(&BIN2BCD8[top*4]));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1044 UBFROMUI(ub+3, UBTOUI(&BIN2BCD8[mid*4]));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1045 UBFROMUI(ub+6, UBTOUI(&BIN2BCD8[rem*4]));
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1046 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1047 else { /* *pa==0 */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1048 UBFROMUI(ub, 0); /* clear 9 BCD8s */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1049 UBFROMUI(ub+4, 0); /* .. */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1050 *(ub+8)=0; /* .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1051 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1052 if (pa==acc) break;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1053 } /* BCD conversion loop */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1054
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1055 num->lsd=ub+8; /* complete the bcdnum .. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1056
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1057 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1058 decShowNum(num, "postmult");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1059 decFloatShow(dfl, "dfl");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1060 decFloatShow(dfr, "dfr");
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1061 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1062 return;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1063 } /* decFiniteMultiply */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1064
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1065 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1066 /* decFloatAbs -- absolute value, heeding NaNs, etc. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1067 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1068 /* result gets the canonicalized df with sign 0 */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1069 /* df is the decFloat to abs */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1070 /* set is the context */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1071 /* returns result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1072 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1073 /* This has the same effect as decFloatPlus unless df is negative, */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1074 /* in which case it has the same effect as decFloatMinus. The */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1075 /* effect is also the same as decFloatCopyAbs except that NaNs are */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1076 /* handled normally (the sign of a NaN is not affected, and an sNaN */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1077 /* will signal) and the result will be canonical. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1078 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1079 decFloat * decFloatAbs(decFloat *result, const decFloat *df,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1080 decContext *set) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1081 if (DFISNAN(df)) return decNaNs(result, df, NULL, set);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1082 decCanonical(result, df); /* copy and check */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1083 DFBYTE(result, 0)&=~0x80; /* zero sign bit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1084 return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1085 } /* decFloatAbs */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1086
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1087 /* ------------------------------------------------------------------ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1088 /* decFloatAdd -- add two decFloats */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1089 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1090 /* result gets the result of adding dfl and dfr: */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1091 /* dfl is the first decFloat (lhs) */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1092 /* dfr is the second decFloat (rhs) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1093 /* set is the context */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1094 /* returns result */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1095 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1096 /* ------------------------------------------------------------------ */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1097 #if QUAD
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1098 /* Table for testing MSDs for fastpath elimination; returns the MSD of */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1099 /* a decDouble or decQuad (top 6 bits tested) ignoring the sign. */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1100 /* Infinities return -32 and NaNs return -128 so that summing the two */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1101 /* MSDs also allows rapid tests for the Specials (see code below). */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1102 const Int DECTESTMSD[64]={
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1103 0, 1, 2, 3, 4, 5, 6, 7, 0, 1, 2, 3, 4, 5, 6, 7,
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1104 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 9, 8, 9, -32, -128,
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1105 0, 1, 2, 3, 4, 5, 6, 7, 0, 1, 2, 3, 4, 5, 6, 7,
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1106 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 9, 8, 9, -32, -128};
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1107 #else
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1108 /* The table for testing MSDs is shared between the modules */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1109 extern const Int DECTESTMSD[64];
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1110 #endif
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1111
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1112 decFloat * decFloatAdd(decFloat *result,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1113 const decFloat *dfl, const decFloat *dfr,
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1114 decContext *set) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1115 bcdnum num; /* for final conversion */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1116 Int bexpl, bexpr; /* left and right biased exponents */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1117 uByte *ub, *us, *ut; /* work */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1118 uInt uiwork; /* for macros */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1119 #if QUAD
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1120 uShort uswork; /* .. */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1121 #endif
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1122
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1123 uInt sourhil, sourhir; /* top words from source decFloats */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1124 /* [valid only through end of */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1125 /* fastpath code -- before swap] */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1126 uInt diffsign; /* non-zero if signs differ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1127 uInt carry; /* carry: 0 or 1 before add loop */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1128 Int overlap; /* coefficient overlap (if full) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1129 Int summ; /* sum of the MSDs */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1130 /* the following buffers hold coefficients with various alignments */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1131 /* (see commentary and diagrams below) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1132 uByte acc[4+2+DECPMAX*3+8];
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1133 uByte buf[4+2+DECPMAX*2];
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1134 uByte *umsd, *ulsd; /* local MSD and LSD pointers */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1135
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1136 #if DECLITEND
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1137 #define CARRYPAT 0x01000000 /* carry=1 pattern */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1138 #else
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1139 #define CARRYPAT 0x00000001 /* carry=1 pattern */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1140 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1141
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1142 /* Start decoding the arguments */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1143 /* The initial exponents are placed into the opposite Ints to */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1144 /* that which might be expected; there are two sets of data to */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1145 /* keep track of (each decFloat and the corresponding exponent), */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1146 /* and this scheme means that at the swap point (after comparing */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1147 /* exponents) only one pair of words needs to be swapped */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1148 /* whichever path is taken (thereby minimising worst-case path). */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1149 /* The calculated exponents will be nonsense when the arguments are */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1150 /* Special, but are not used in that path */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1151 sourhil=DFWORD(dfl, 0); /* LHS top word */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1152 summ=DECTESTMSD[sourhil>>26]; /* get first MSD for testing */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1153 bexpr=DECCOMBEXP[sourhil>>26]; /* get exponent high bits (in place) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1154 bexpr+=GETECON(dfl); /* .. + continuation */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1155
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1156 sourhir=DFWORD(dfr, 0); /* RHS top word */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1157 summ+=DECTESTMSD[sourhir>>26]; /* sum MSDs for testing */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1158 bexpl=DECCOMBEXP[sourhir>>26];
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1159 bexpl+=GETECON(dfr);
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1160
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1161 /* here bexpr has biased exponent from lhs, and vice versa */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1162
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1163 diffsign=(sourhil^sourhir)&DECFLOAT_Sign;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1164
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1165 /* now determine whether to take a fast path or the full-function */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1166 /* slow path. The slow path must be taken when: */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1167 /* -- both numbers are finite, and: */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1168 /* the exponents are different, or */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1169 /* the signs are different, or */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1170 /* the sum of the MSDs is >8 (hence might overflow) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1171 /* specialness and the sum of the MSDs can be tested at once using */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1172 /* the summ value just calculated, so the test for specials is no */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1173 /* longer on the worst-case path (as of 3.60) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1174
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1175 if (summ<=8) { /* MSD+MSD is good, or there is a special */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1176 if (summ<0) { /* there is a special */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1177 /* Inf+Inf would give -64; Inf+finite is -32 or higher */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1178 if (summ<-64) return decNaNs(result, dfl, dfr, set); /* one or two NaNs */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1179 /* two infinities with different signs is invalid */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1180 if (summ==-64 && diffsign) return decInvalid(result, set);
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1181 if (DFISINF(dfl)) return decInfinity(result, dfl); /* LHS is infinite */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1182 return decInfinity(result, dfr); /* RHS must be Inf */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1183 }
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1184 /* Here when both arguments are finite; fast path is possible */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1185 /* (currently only for aligned and same-sign) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1186 if (bexpr==bexpl && !diffsign) {
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1187 uInt tac[DECLETS+1]; /* base-1000 coefficient */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1188 uInt encode; /* work */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1189
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1190 /* Get one coefficient as base-1000 and add the other */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1191 GETCOEFFTHOU(dfl, tac); /* least-significant goes to [0] */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1192 ADDCOEFFTHOU(dfr, tac);
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1193 /* here the sum of the MSDs (plus any carry) will be <10 due to */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1194 /* the fastpath test earlier */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1195
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1196 /* construct the result; low word is the same for both formats */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1197 encode =BIN2DPD[tac[0]];
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1198 encode|=BIN2DPD[tac[1]]<<10;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1199 encode|=BIN2DPD[tac[2]]<<20;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1200 encode|=BIN2DPD[tac[3]]<<30;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1201 DFWORD(result, (DECBYTES/4)-1)=encode;
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1202
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1203 /* collect next two declets (all that remains, for Double) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1204 encode =BIN2DPD[tac[3]]>>2;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1205 encode|=BIN2DPD[tac[4]]<<8;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1206
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1207 #if QUAD
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1208 /* complete and lay out middling words */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1209 encode|=BIN2DPD[tac[5]]<<18;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1210 encode|=BIN2DPD[tac[6]]<<28;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1211 DFWORD(result, 2)=encode;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1212
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1213 encode =BIN2DPD[tac[6]]>>4;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1214 encode|=BIN2DPD[tac[7]]<<6;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1215 encode|=BIN2DPD[tac[8]]<<16;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1216 encode|=BIN2DPD[tac[9]]<<26;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1217 DFWORD(result, 1)=encode;
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1218
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1219 /* and final two declets */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1220 encode =BIN2DPD[tac[9]]>>6;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1221 encode|=BIN2DPD[tac[10]]<<4;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1222 #endif
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1223
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1224 /* add exponent continuation and sign (from either argument) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1225 encode|=sourhil & (ECONMASK | DECFLOAT_Sign);
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1226
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1227 /* create lookup index = MSD + top two bits of biased exponent <<4 */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1228 tac[DECLETS]|=(bexpl>>DECECONL)<<4;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1229 encode|=DECCOMBFROM[tac[DECLETS]]; /* add constructed combination field */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1230 DFWORD(result, 0)=encode; /* complete */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1231
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1232 /* decFloatShow(result, ">"); */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1233 return result;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1234 } /* fast path OK */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1235 /* drop through to slow path */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1236 } /* low sum or Special(s) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1237
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1238 /* Slow path required -- arguments are finite and might overflow, */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1239 /* or require alignment, or might have different signs */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1240
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1241 /* now swap either exponents or argument pointers */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1242 if (bexpl<=bexpr) {
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1243 /* original left is bigger */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1244 Int bexpswap=bexpl;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1245 bexpl=bexpr;
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1246 bexpr=bexpswap;
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1247 /* printf("left bigger\n"); */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1248 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1249 else {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1250 const decFloat *dfswap=dfl;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1251 dfl=dfr;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1252 dfr=dfswap;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1253 /* printf("right bigger\n"); */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1254 }
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1255 /* [here dfl and bexpl refer to the datum with the larger exponent, */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1256 /* of if the exponents are equal then the original LHS argument] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1257
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1258 /* if lhs is zero then result will be the rhs (now known to have */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1259 /* the smaller exponent), which also may need to be tested for zero */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1260 /* for the weird IEEE 754 sign rules */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1261 if (DFISZERO(dfl)) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1262 decCanonical(result, dfr); /* clean copy */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1263 /* "When the sum of two operands with opposite signs is */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1264 /* exactly zero, the sign of that sum shall be '+' in all */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1265 /* rounding modes except round toward -Infinity, in which */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1266 /* mode that sign shall be '-'." */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1267 if (diffsign && DFISZERO(result)) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1268 DFWORD(result, 0)&=~DECFLOAT_Sign; /* assume sign 0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1269 if (set->round==DEC_ROUND_FLOOR) DFWORD(result, 0)|=DECFLOAT_Sign;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1270 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1271 return result;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1272 } /* numfl is zero */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1273 /* [here, LHS is non-zero; code below assumes that] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1274
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1275 /* Coefficients layout during the calculations to follow: */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1276 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1277 /* Overlap case: */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1278 /* +------------------------------------------------+ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1279 /* acc: |0000| coeffa | tail B | | */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1280 /* +------------------------------------------------+ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1281 /* buf: |0000| pad0s | coeffb | | */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1282 /* +------------------------------------------------+ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1283 /* */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1284 /* Touching coefficients or gap: */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1285 /* +------------------------------------------------+ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1286 /* acc: |0000| coeffa | gap | coeffb | */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1287 /* +------------------------------------------------+ */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1288 /* [buf not used or needed; gap clamped to Pmax] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1289
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1290 /* lay out lhs coefficient into accumulator; this starts at acc+4 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1291 /* for decDouble or acc+6 for decQuad so the LSD is word- */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1292 /* aligned; the top word gap is there only in case a carry digit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1293 /* is prefixed after the add -- it does not need to be zeroed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1294 #if DOUBLE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1295 #define COFF 4 /* offset into acc */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1296 #elif QUAD
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1297 UBFROMUS(acc+4, 0); /* prefix 00 */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1298 #define COFF 6 /* offset into acc */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1299 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1300
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1301 GETCOEFF(dfl, acc+COFF); /* decode from decFloat */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1302 ulsd=acc+COFF+DECPMAX-1;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1303 umsd=acc+4; /* [having this here avoids */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1304
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1305 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1306 {bcdnum tum;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1307 tum.msd=umsd;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1308 tum.lsd=ulsd;
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1309 tum.exponent=bexpl-DECBIAS;
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1310 tum.sign=DFWORD(dfl, 0) & DECFLOAT_Sign;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1311 decShowNum(&tum, "dflx");}
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1312 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1313
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1314 /* if signs differ, take ten's complement of lhs (here the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1315 /* coefficient is subtracted from all-nines; the 1 is added during */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1316 /* the later add cycle -- zeros to the right do not matter because */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1317 /* the complement of zero is zero); these are fixed-length inverts */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1318 /* where the lsd is known to be at a 4-byte boundary (so no borrow */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1319 /* possible) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1320 carry=0; /* assume no carry */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1321 if (diffsign) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1322 carry=CARRYPAT; /* for +1 during add */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1323 UBFROMUI(acc+ 4, 0x09090909-UBTOUI(acc+ 4));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1324 UBFROMUI(acc+ 8, 0x09090909-UBTOUI(acc+ 8));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1325 UBFROMUI(acc+12, 0x09090909-UBTOUI(acc+12));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1326 UBFROMUI(acc+16, 0x09090909-UBTOUI(acc+16));
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1327 #if QUAD
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1328 UBFROMUI(acc+20, 0x09090909-UBTOUI(acc+20));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1329 UBFROMUI(acc+24, 0x09090909-UBTOUI(acc+24));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1330 UBFROMUI(acc+28, 0x09090909-UBTOUI(acc+28));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1331 UBFROMUI(acc+32, 0x09090909-UBTOUI(acc+32));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1332 UBFROMUI(acc+36, 0x09090909-UBTOUI(acc+36));
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1333 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1334 } /* diffsign */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1335
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1336 /* now process the rhs coefficient; if it cannot overlap lhs then */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1337 /* it can be put straight into acc (with an appropriate gap, if */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1338 /* needed) because no actual addition will be needed (except */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1339 /* possibly to complete ten's complement) */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1340 overlap=DECPMAX-(bexpl-bexpr);
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1341 #if DECTRACE
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1342 printf("exps: %ld %ld\n", (LI)(bexpl-DECBIAS), (LI)(bexpr-DECBIAS));
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1343 printf("Overlap=%ld carry=%08lx\n", (LI)overlap, (LI)carry);
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1344 #endif
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1345
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1346 if (overlap<=0) { /* no overlap possible */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1347 uInt gap; /* local work */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1348 /* since a full addition is not needed, a ten's complement */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1349 /* calculation started above may need to be completed */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1350 if (carry) {
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1351 for (ub=ulsd; *ub==9; ub--) *ub=0;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1352 *ub+=1;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1353 carry=0; /* taken care of */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1354 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1355 /* up to DECPMAX-1 digits of the final result can extend down */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1356 /* below the LSD of the lhs, so if the gap is >DECPMAX then the */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1357 /* rhs will be simply sticky bits. In this case the gap is */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1358 /* clamped to DECPMAX and the exponent adjusted to suit [this is */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1359 /* safe because the lhs is non-zero]. */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1360 gap=-overlap;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1361 if (gap>DECPMAX) {
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1362 bexpr+=gap-1;
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1363 gap=DECPMAX;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1364 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1365 ub=ulsd+gap+1; /* where MSD will go */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1366 /* Fill the gap with 0s; note that there is no addition to do */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1367 ut=acc+COFF+DECPMAX; /* start of gap */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1368 for (; ut<ub; ut+=4) UBFROMUI(ut, 0); /* mind the gap */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1369 if (overlap<-DECPMAX) { /* gap was > DECPMAX */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1370 *ub=(uByte)(!DFISZERO(dfr)); /* make sticky digit */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1371 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1372 else { /* need full coefficient */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1373 GETCOEFF(dfr, ub); /* decode from decFloat */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1374 ub+=DECPMAX-1; /* new LSD... */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1375 }
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1376 ulsd=ub; /* save new LSD */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1377 } /* no overlap possible */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1378
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1379 else { /* overlap>0 */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1380 /* coefficients overlap (perhaps completely, although also */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1381 /* perhaps only where zeros) */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1382 if (overlap==DECPMAX) { /* aligned */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1383 ub=buf+COFF; /* where msd will go */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1384 #if QUAD
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1385 UBFROMUS(buf+4, 0); /* clear quad's 00 */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1386 #endif
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1387 GETCOEFF(dfr, ub); /* decode from decFloat */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1388 }
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1389 else { /* unaligned */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1390 ub=buf+COFF+DECPMAX-overlap; /* where MSD will go */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1391 /* Fill the prefix gap with 0s; 8 will cover most common */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1392 /* unalignments, so start with direct assignments (a loop is */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1393 /* then used for any remaining -- the loop (and the one in a */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1394 /* moment) is not then on the critical path because the number */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1395 /* of additions is reduced by (at least) two in this case) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1396 UBFROMUI(buf+4, 0); /* [clears decQuad 00 too] */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1397 UBFROMUI(buf+8, 0);
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1398 if (ub>buf+12) {
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1399 ut=buf+12; /* start any remaining */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1400 for (; ut<ub; ut+=4) UBFROMUI(ut, 0); /* fill them */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1401 }
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1402 GETCOEFF(dfr, ub); /* decode from decFloat */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1403
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1404 /* now move tail of rhs across to main acc; again use direct */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1405 /* copies for 8 digits-worth */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1406 UBFROMUI(acc+COFF+DECPMAX, UBTOUI(buf+COFF+DECPMAX));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1407 UBFROMUI(acc+COFF+DECPMAX+4, UBTOUI(buf+COFF+DECPMAX+4));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1408 if (buf+COFF+DECPMAX+8<ub+DECPMAX) {
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1409 us=buf+COFF+DECPMAX+8; /* source */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1410 ut=acc+COFF+DECPMAX+8; /* target */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1411 for (; us<ub+DECPMAX; us+=4, ut+=4) UBFROMUI(ut, UBTOUI(us));
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1412 }
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1413 } /* unaligned */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1414
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1415 ulsd=acc+(ub-buf+DECPMAX-1); /* update LSD pointer */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1416
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1417 /* Now do the add of the non-tail; this is all nicely aligned, */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1418 /* and is over a multiple of four digits (because for Quad two */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1419 /* zero digits were added on the left); words in both acc and */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1420 /* buf (buf especially) will often be zero */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1421 /* [byte-by-byte add, here, is about 15% slower total effect than */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1422 /* the by-fours] */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1423
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1424 /* Now effect the add; this is harder on a little-endian */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1425 /* machine as the inter-digit carry cannot use the usual BCD */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1426 /* addition trick because the bytes are loaded in the wrong order */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1427 /* [this loop could be unrolled, but probably scarcely worth it] */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1428
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1429 ut=acc+COFF+DECPMAX-4; /* target LSW (acc) */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1430 us=buf+COFF+DECPMAX-4; /* source LSW (buf, to add to acc) */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1431
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1432 #if !DECLITEND
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1433 for (; ut>=acc+4; ut-=4, us-=4) { /* big-endian add loop */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1434 /* bcd8 add */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1435 carry+=UBTOUI(us); /* rhs + carry */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1436 if (carry==0) continue; /* no-op */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1437 carry+=UBTOUI(ut); /* lhs */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1438 /* Big-endian BCD adjust (uses internal carry) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1439 carry+=0x76f6f6f6; /* note top nibble not all bits */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1440 /* apply BCD adjust and save */
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1441 UBFROMUI(ut, (carry & 0x0f0f0f0f) - ((carry & 0x60606060)>>4));
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1442 carry>>=31; /* true carry was at far left */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1443 } /* add loop */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1444 #else
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1445 for (; ut>=acc+4; ut-=4, us-=4) { /* little-endian add loop */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1446 /* bcd8 add */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1447 carry+=UBTOUI(us); /* rhs + carry */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1448 if (carry==0) continue; /* no-op [common if unaligned] */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1449 carry+=UBTOUI(ut); /* lhs */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1450 /* Little-endian BCD adjust; inter-digit carry must be manual */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1451 /* because the lsb from the array will be in the most-significant */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1452 /* byte of carry */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1453 carry+=0x76767676; /* note no inter-byte carries */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1454 carry+=(carry & 0x80000000)>>15;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1455 carry+=(carry & 0x00800000)>>15;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1456 carry+=(carry & 0x00008000)>>15;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1457 carry-=(carry & 0x60606060)>>4; /* BCD adjust back */
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1458 UBFROMUI(ut, carry & 0x0f0f0f0f); /* clear debris and save */
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1459 /* here, final carry-out bit is at 0x00000080; move it ready */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1460 /* for next word-add (i.e., to 0x01000000) */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1461 carry=(carry & 0x00000080)<<17;
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1462 } /* add loop */
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1463 #endif
55
77e2b8dfacca update it from 4.4.3 to 4.5.0
ryoma <e075725@ie.u-ryukyu.ac.jp>
parents: 0
diff changeset
1464
0
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1465 #if DECTRACE
a06113de4d67 first commit
kent <kent@cr.ie.u-ryukyu.ac.jp>
parents:
diff changeset
1466 {bcdnum tum;
a06113de4d67 first commit