merge from gcc
[platform/upstream/binutils.git] / libdecnumber / decCommon.c
1 /* Common code for fixed-size types in the decNumber C Library.
2    Copyright (C) 2007 Free Software Foundation, Inc.
3    Contributed by IBM Corporation.  Author Mike Cowlishaw.
4
5    This file is part of GCC.
6
7    GCC is free software; you can redistribute it and/or modify it under
8    the terms of the GNU General Public License as published by the Free
9    Software Foundation; either version 2, or (at your option) any later
10    version.
11
12    In addition to the permissions in the GNU General Public License,
13    the Free Software Foundation gives you unlimited permission to link
14    the compiled version of this file into combinations with other
15    programs, and to distribute those combinations without any
16    restriction coming from the use of this file.  (The General Public
17    License restrictions do apply in other respects; for example, they
18    cover modification of the file, and distribution when not linked
19    into a combine executable.)
20
21    GCC is distributed in the hope that it will be useful, but WITHOUT ANY
22    WARRANTY; without even the implied warranty of MERCHANTABILITY or
23    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
24    for more details.
25
26    You should have received a copy of the GNU General Public License
27    along with GCC; see the file COPYING.  If not, write to the Free
28    Software Foundation, 51 Franklin Street, Fifth Floor, Boston, MA
29    02110-1301, USA.  */
30
31 /* ------------------------------------------------------------------ */
32 /* decCommon.c -- common code for all three fixed-size types          */
33 /* ------------------------------------------------------------------ */
34 /* This module comprises code that is shared between all the formats  */
35 /* (decSingle, decDouble, and decQuad); it includes set and extract   */
36 /* of format components, widening, narrowing, and string conversions. */
37 /*                                                                    */
38 /* Unlike decNumber, parameterization takes place at compile time     */
39 /* rather than at runtime.  The parameters are set in the decDouble.c */
40 /* (etc.) files, which then include this one to produce the compiled  */
41 /* code.  The functions here, therefore, are code shared between      */
42 /* multiple formats.                                                  */
43 /* ------------------------------------------------------------------ */
44 /* Names here refer to decFloat rather than to decDouble, etc., and */
45 /* the functions are in strict alphabetical order. */
46 /* Constants, tables, and debug function(s) are included only for QUAD */
47 /* (which will always be compiled if DOUBLE or SINGLE are used). */
48 /* */
49 /* Whenever a decContext is used, only the status may be set (using */
50 /* OR) or the rounding mode read; all other fields are ignored and */
51 /* untouched. */
52
53 #include "decCommonSymbols.h"
54
55 /* names for simpler testing and default context */
56 #if DECPMAX==7
57   #define SINGLE     1
58   #define DOUBLE     0
59   #define QUAD       0
60   #define DEFCONTEXT DEC_INIT_DECIMAL32
61 #elif DECPMAX==16
62   #define SINGLE     0
63   #define DOUBLE     1
64   #define QUAD       0
65   #define DEFCONTEXT DEC_INIT_DECIMAL64
66 #elif DECPMAX==34
67   #define SINGLE     0
68   #define DOUBLE     0
69   #define QUAD       1
70   #define DEFCONTEXT DEC_INIT_DECIMAL128
71 #else
72   #error Unexpected DECPMAX value
73 #endif
74
75 /* Assertions */
76
77 #if DECPMAX!=7 && DECPMAX!=16 && DECPMAX!=34
78   #error Unexpected Pmax (DECPMAX) value for this module
79 #endif
80
81 /* Assert facts about digit characters, etc. */
82 #if ('9'&0x0f)!=9
83   #error This module assumes characters are of the form 0b....nnnn
84   /* where .... are don't care 4 bits and nnnn is 0000 through 1001 */
85 #endif
86 #if ('9'&0xf0)==('.'&0xf0)
87   #error This module assumes '.' has a different mask than a digit
88 #endif
89
90 /* Assert ToString lay-out conditions */
91 #if DECSTRING<DECPMAX+9
92   #error ToString needs at least 8 characters for lead-in and dot
93 #endif
94 #if DECPMAX+DECEMAXD+5 > DECSTRING
95   #error Exponent form can be too long for ToString to lay out safely
96 #endif
97 #if DECEMAXD > 4
98   #error Exponent form is too long for ToString to lay out
99   /* Note: code for up to 9 digits exists in archives [decOct] */
100 #endif
101
102 /* Private functions used here and possibly in decBasic.c, etc. */
103 static decFloat * decFinalize(decFloat *, bcdnum *, decContext *);
104 static Flag decBiStr(const char *, const char *, const char *);
105
106 /* Macros and private tables; those which are not format-dependent    */
107 /* are only included if decQuad is being built.                       */
108
109 /* ------------------------------------------------------------------ */
110 /* Combination field lookup tables (uInts to save measurable work)    */
111 /*                                                                    */
112 /*   DECCOMBEXP  - 2 most-significant-bits of exponent (00, 01, or    */
113 /*                 10), shifted left for format, or DECFLOAT_Inf/NaN  */
114 /*   DECCOMBWEXP - The same, for the next-wider format (unless QUAD)  */
115 /*   DECCOMBMSD  - 4-bit most-significant-digit                       */
116 /*                 [0 if the index is a special (Infinity or NaN)]    */
117 /*   DECCOMBFROM - 5-bit combination field from EXP top bits and MSD  */
118 /*                 (placed in uInt so no shift is needed)             */
119 /*                                                                    */
120 /* DECCOMBEXP, DECCOMBWEXP, and DECCOMBMSD are indexed by the sign    */
121 /*   and 5-bit combination field (0-63, the second half of the table  */
122 /*   identical to the first half)                                     */
123 /* DECCOMBFROM is indexed by expTopTwoBits*16 + msd                   */
124 /*                                                                    */
125 /* DECCOMBMSD and DECCOMBFROM are not format-dependent and so are     */
126 /* only included once, when QUAD is being built                       */
127 /* ------------------------------------------------------------------ */
128 static const uInt DECCOMBEXP[64]={
129   0, 0, 0, 0, 0, 0, 0, 0,
130   1<<DECECONL, 1<<DECECONL, 1<<DECECONL, 1<<DECECONL,
131   1<<DECECONL, 1<<DECECONL, 1<<DECECONL, 1<<DECECONL,
132   2<<DECECONL, 2<<DECECONL, 2<<DECECONL, 2<<DECECONL,
133   2<<DECECONL, 2<<DECECONL, 2<<DECECONL, 2<<DECECONL,
134   0,           0,           1<<DECECONL, 1<<DECECONL,
135   2<<DECECONL, 2<<DECECONL, DECFLOAT_Inf, DECFLOAT_NaN,
136   0, 0, 0, 0, 0, 0, 0, 0,
137   1<<DECECONL, 1<<DECECONL, 1<<DECECONL, 1<<DECECONL,
138   1<<DECECONL, 1<<DECECONL, 1<<DECECONL, 1<<DECECONL,
139   2<<DECECONL, 2<<DECECONL, 2<<DECECONL, 2<<DECECONL,
140   2<<DECECONL, 2<<DECECONL, 2<<DECECONL, 2<<DECECONL,
141   0,           0,           1<<DECECONL, 1<<DECECONL,
142   2<<DECECONL, 2<<DECECONL, DECFLOAT_Inf, DECFLOAT_NaN};
143 #if !QUAD
144 static const uInt DECCOMBWEXP[64]={
145   0, 0, 0, 0, 0, 0, 0, 0,
146   1<<DECWECONL, 1<<DECWECONL, 1<<DECWECONL, 1<<DECWECONL,
147   1<<DECWECONL, 1<<DECWECONL, 1<<DECWECONL, 1<<DECWECONL,
148   2<<DECWECONL, 2<<DECWECONL, 2<<DECWECONL, 2<<DECWECONL,
149   2<<DECWECONL, 2<<DECWECONL, 2<<DECWECONL, 2<<DECWECONL,
150   0,            0,            1<<DECWECONL, 1<<DECWECONL,
151   2<<DECWECONL, 2<<DECWECONL, DECFLOAT_Inf, DECFLOAT_NaN,
152   0, 0, 0, 0, 0, 0, 0, 0,
153   1<<DECWECONL, 1<<DECWECONL, 1<<DECWECONL, 1<<DECWECONL,
154   1<<DECWECONL, 1<<DECWECONL, 1<<DECWECONL, 1<<DECWECONL,
155   2<<DECWECONL, 2<<DECWECONL, 2<<DECWECONL, 2<<DECWECONL,
156   2<<DECWECONL, 2<<DECWECONL, 2<<DECWECONL, 2<<DECWECONL,
157   0,            0,            1<<DECWECONL, 1<<DECWECONL,
158   2<<DECWECONL, 2<<DECWECONL, DECFLOAT_Inf, DECFLOAT_NaN};
159 #endif
160
161 #if QUAD
162 const uInt DECCOMBMSD[64]={
163   0, 1, 2, 3, 4, 5, 6, 7, 0, 1, 2, 3, 4, 5, 6, 7,
164   0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 9, 8, 9, 0, 1,
165   0, 1, 2, 3, 4, 5, 6, 7, 0, 1, 2, 3, 4, 5, 6, 7,
166   0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 9, 8, 9, 0, 0};
167
168 const uInt DECCOMBFROM[48]={
169   0x00000000, 0x04000000, 0x08000000, 0x0C000000, 0x10000000, 0x14000000,
170   0x18000000, 0x1C000000, 0x60000000, 0x64000000, 0x00000000, 0x00000000,
171   0x00000000, 0x00000000, 0x00000000, 0x00000000, 0x20000000, 0x24000000,
172   0x28000000, 0x2C000000, 0x30000000, 0x34000000, 0x38000000, 0x3C000000,
173   0x68000000, 0x6C000000, 0x00000000, 0x00000000, 0x00000000, 0x00000000,
174   0x00000000, 0x00000000, 0x40000000, 0x44000000, 0x48000000, 0x4C000000,
175   0x50000000, 0x54000000, 0x58000000, 0x5C000000, 0x70000000, 0x74000000,
176   0x00000000, 0x00000000, 0x00000000, 0x00000000, 0x00000000, 0x00000000};
177
178 /* ------------------------------------------------------------------ */
179 /* Request and include the tables to use for conversions              */
180 /* ------------------------------------------------------------------ */
181 #define DEC_BCD2DPD  1        /* 0-0x999 -> DPD */
182 #define DEC_BIN2DPD  1        /* 0-999 -> DPD */
183 #define DEC_BIN2BCD8 1        /* 0-999 -> ddd, len */
184 #define DEC_DPD2BCD8 1        /* DPD -> ddd, len */
185 #define DEC_DPD2BIN  1        /* DPD -> 0-999 */
186 #define DEC_DPD2BINK 1        /* DPD -> 0-999000 */
187 #define DEC_DPD2BINM 1        /* DPD -> 0-999000000 */
188 #include "decDPD.h"           /* source of the lookup tables */
189
190 #endif
191
192 /* ----------------------------------------------------------------- */
193 /* decBiStr -- compare string with pairwise options                  */
194 /*                                                                   */
195 /*   targ is the string to compare                                   */
196 /*   str1 is one of the strings to compare against (length may be 0) */
197 /*   str2 is the other; it must be the same length as str1           */
198 /*                                                                   */
199 /*   returns 1 if strings compare equal, (that is, targ is the same  */
200 /*   length as str1 and str2, and each character of targ is in one   */
201 /*   of str1 or str2 in the corresponding position), or 0 otherwise  */
202 /*                                                                   */
203 /* This is used for generic caseless compare, including the awkward  */
204 /* case of the Turkish dotted and dotless Is.  Use as (for example): */
205 /*   if (decBiStr(test, "mike", "MIKE")) ...                         */
206 /* ----------------------------------------------------------------- */
207 static Flag decBiStr(const char *targ, const char *str1, const char *str2) {
208   for (;;targ++, str1++, str2++) {
209     if (*targ!=*str1 && *targ!=*str2) return 0;
210     /* *targ has a match in one (or both, if terminator) */
211     if (*targ=='\0') break;
212     } /* forever */
213   return 1;
214   } /* decBiStr */
215
216 /* ------------------------------------------------------------------ */
217 /* decFinalize -- adjust and store a final result                     */
218 /*                                                                    */
219 /*  df  is the decFloat format number which gets the final result     */
220 /*  num is the descriptor of the number to be checked and encoded     */
221 /*         [its values, including the coefficient, may be modified]   */
222 /*  set is the context to use                                         */
223 /*  returns df                                                        */
224 /*                                                                    */
225 /* The num descriptor may point to a bcd8 string of any length; this  */
226 /* string may have leading insignificant zeros.  If it has more than  */
227 /* DECPMAX digits then the final digit can be a round-for-reround     */
228 /* digit (i.e., it may include a sticky bit residue).                 */
229 /*                                                                    */
230 /* The exponent (q) may be one of the codes for a special value and   */
231 /* can be up to 999999999 for conversion from string.                 */
232 /*                                                                    */
233 /* No error is possible, but Inexact, Underflow, and/or Overflow may  */
234 /* be set.                                                            */
235 /* ------------------------------------------------------------------ */
236 /* Constant whose size varies with format; also the check for surprises */
237 static uByte allnines[DECPMAX]=
238 #if SINGLE
239   {9, 9, 9, 9, 9, 9, 9};
240 #elif DOUBLE
241   {9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9};
242 #elif QUAD
243   {9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9,
244    9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9, 9};
245 #endif
246
247 static decFloat * decFinalize(decFloat *df, bcdnum *num,
248                               decContext *set) {
249   uByte *ub;                  /* work */
250   uInt   dpd;                 /* .. */
251   uByte *umsd=num->msd;       /* local copy */
252   uByte *ulsd=num->lsd;       /* .. */
253   uInt   encode;              /* encoding accumulator */
254   Int    length;              /* coefficient length */
255
256   #if DECCHECK
257   Int clen=ulsd-umsd+1;
258   #if QUAD
259     #define COEXTRA 2                        /* extra-long coefficent */
260   #else
261     #define COEXTRA 0
262   #endif
263   if (clen<1 || clen>DECPMAX*3+2+COEXTRA)
264     printf("decFinalize: suspect coefficient [length=%ld]\n", (LI)clen);
265   if (num->sign!=0 && num->sign!=DECFLOAT_Sign)
266     printf("decFinalize: bad sign [%08lx]\n", (LI)num->sign);
267   if (!EXPISSPECIAL(num->exponent)
268       && (num->exponent>1999999999 || num->exponent<-1999999999))
269     printf("decFinalize: improbable exponent [%ld]\n", (LI)num->exponent);
270   /* decShowNum(num, "final"); */
271   #endif
272
273   /* A special will have an 'exponent' which is very positive and a */
274   /* coefficient < DECPMAX */
275   length=(uInt)(ulsd-umsd+1);                /* coefficient length */
276
277   if (!NUMISSPECIAL(num)) {
278     Int   drop;                              /* digits to be dropped */
279     /* skip leading insignificant zeros to calculate an exact length */
280     /* [this is quite expensive] */
281     if (*umsd==0) {
282       for (; UINTAT(umsd)==0 && umsd+3<ulsd;) umsd+=4;
283       for (; *umsd==0 && umsd<ulsd;) umsd++;
284       length=ulsd-umsd+1;                    /* recalculate */
285       }
286     drop=MAXI(length-DECPMAX, DECQTINY-num->exponent);
287     /* drop can now be > digits for bottom-clamp (subnormal) cases */
288     if (drop>0) {                            /* rounding needed */
289       /* (decFloatQuantize has very similar code to this, so any */
290       /* changes may need to be made there, too) */
291       uByte *roundat;                        /* -> re-round digit */
292       uByte reround;                         /* reround value */
293       /* printf("Rounding; drop=%ld\n", (LI)drop); */
294
295       num->exponent+=drop;                   /* always update exponent */
296
297       /* Three cases here: */
298       /*   1. new LSD is in coefficient (almost always) */
299       /*   2. new LSD is digit to left of coefficient (so MSD is */
300       /*      round-for-reround digit) */
301       /*   3. new LSD is to left of case 2 (whole coefficient is sticky) */
302       /* [duplicate check-stickies code to save a test] */
303       /* [by-digit check for stickies as runs of zeros are rare] */
304       if (drop<length) {                     /* NB lengths not addresses */
305         roundat=umsd+length-drop;
306         reround=*roundat;
307         for (ub=roundat+1; ub<=ulsd; ub++) {
308           if (*ub!=0) {                      /* non-zero to be discarded */
309             reround=DECSTICKYTAB[reround];   /* apply sticky bit */
310             break;                           /* [remainder don't-care] */
311             }
312           } /* check stickies */
313         ulsd=roundat-1;                      /* new LSD */
314         }
315        else {                                /* edge case */
316         if (drop==length) {
317           roundat=umsd;
318           reround=*roundat;
319           }
320          else {
321           roundat=umsd-1;
322           reround=0;
323           }
324         for (ub=roundat+1; ub<=ulsd; ub++) {
325           if (*ub!=0) {                      /* non-zero to be discarded */
326             reround=DECSTICKYTAB[reround];   /* apply sticky bit */
327             break;                           /* [remainder don't-care] */
328             }
329           } /* check stickies */
330         *umsd=0;                             /* coefficient is a 0 */
331         ulsd=umsd;                           /* .. */
332         }
333
334       if (reround!=0) {                      /* discarding non-zero */
335         uInt bump=0;
336         set->status|=DEC_Inexact;
337         /* if adjusted exponent [exp+digits-1] is < EMIN then num is */
338         /* subnormal -- so raise Underflow */
339         if (num->exponent<DECEMIN && (num->exponent+(ulsd-umsd+1)-1)<DECEMIN)
340           set->status|=DEC_Underflow;
341
342         /* next decide whether increment of the coefficient is needed */
343         if (set->round==DEC_ROUND_HALF_EVEN) {    /* fastpath slowest case */
344           if (reround>5) bump=1;                  /* >0.5 goes up */
345            else if (reround==5)                   /* exactly 0.5000 .. */
346             bump=*ulsd & 0x01;                    /* .. up iff [new] lsd is odd */
347           } /* r-h-e */
348          else switch (set->round) {
349           case DEC_ROUND_DOWN: {
350             /* no change */
351             break;} /* r-d */
352           case DEC_ROUND_HALF_DOWN: {
353             if (reround>5) bump=1;
354             break;} /* r-h-d */
355           case DEC_ROUND_HALF_UP: {
356             if (reround>=5) bump=1;
357             break;} /* r-h-u */
358           case DEC_ROUND_UP: {
359             if (reround>0) bump=1;
360             break;} /* r-u */
361           case DEC_ROUND_CEILING: {
362             /* same as _UP for positive numbers, and as _DOWN for negatives */
363             if (!num->sign && reround>0) bump=1;
364             break;} /* r-c */
365           case DEC_ROUND_FLOOR: {
366             /* same as _UP for negative numbers, and as _DOWN for positive */
367             /* [negative reround cannot occur on 0] */
368             if (num->sign && reround>0) bump=1;
369             break;} /* r-f */
370           case DEC_ROUND_05UP: {
371             if (reround>0) { /* anything out there is 'sticky' */
372               /* bump iff lsd=0 or 5; this cannot carry so it could be */
373               /* effected immediately with no bump -- but the code */
374               /* is clearer if this is done the same way as the others */
375               if (*ulsd==0 || *ulsd==5) bump=1;
376               }
377             break;} /* r-r */
378           default: {      /* e.g., DEC_ROUND_MAX */
379             set->status|=DEC_Invalid_context;
380             #if DECCHECK
381             printf("Unknown rounding mode: %ld\n", (LI)set->round);
382             #endif
383             break;}
384           } /* switch (not r-h-e) */
385         /* printf("ReRound: %ld  bump: %ld\n", (LI)reround, (LI)bump); */
386
387         if (bump!=0) {                       /* need increment */
388           /* increment the coefficient; this might end up with 1000... */
389           /* (after the all nines case) */
390           ub=ulsd;
391           for(; ub-3>=umsd && UINTAT(ub-3)==0x09090909; ub-=4) UINTAT(ub-3)=0;
392           /* [note ub could now be to left of msd, and it is not safe */
393           /* to write to the the left of the msd] */
394           /* now at most 3 digits left to non-9 (usually just the one) */
395           for (; ub>=umsd; *ub=0, ub--) {
396             if (*ub==9) continue;            /* carry */
397             *ub+=1;
398             break;
399             }
400           if (ub<umsd) {                     /* had all-nines */
401             *umsd=1;                         /* coefficient to 1000... */
402             /* usually the 1000... coefficient can be used as-is */
403             if ((ulsd-umsd+1)==DECPMAX) {
404               num->exponent++;
405               }
406              else {
407               /* if coefficient is shorter than Pmax then num is */
408               /* subnormal, so extend it; this is safe as drop>0 */
409               /* (or, if the coefficient was supplied above, it could */
410               /* not be 9); this may make the result normal. */
411               ulsd++;
412               *ulsd=0;
413               /* [exponent unchanged] */
414               #if DECCHECK
415               if (num->exponent!=DECQTINY) /* sanity check */
416                 printf("decFinalize: bad all-nines extend [^%ld, %ld]\n",
417                        (LI)num->exponent, (LI)(ulsd-umsd+1));
418               #endif
419               } /* subnormal extend */
420             } /* had all-nines */
421           } /* bump needed */
422         } /* inexact rounding */
423
424       length=ulsd-umsd+1;               /* recalculate (may be <DECPMAX) */
425       } /* need round (drop>0) */
426
427     /* The coefficient will now fit and has final length unless overflow */
428     /* decShowNum(num, "rounded"); */
429
430     /* if exponent is >=emax may have to clamp, overflow, or fold-down */
431     if (num->exponent>DECEMAX-(DECPMAX-1)) { /* is edge case */
432       /* printf("overflow checks...\n"); */
433       if (*ulsd==0 && ulsd==umsd) {     /* have zero */
434         num->exponent=DECEMAX-(DECPMAX-1); /* clamp to max */
435         }
436        else if ((num->exponent+length-1)>DECEMAX) { /* > Nmax */
437         /* Overflow -- these could go straight to encoding, here, but */
438         /* instead num is adjusted to keep the code cleaner */
439         Flag needmax=0;                 /* 1 for finite result */
440         set->status|=(DEC_Overflow | DEC_Inexact);
441         switch (set->round) {
442           case DEC_ROUND_DOWN: {
443             needmax=1;                  /* never Infinity */
444             break;} /* r-d */
445           case DEC_ROUND_05UP: {
446             needmax=1;                  /* never Infinity */
447             break;} /* r-05 */
448           case DEC_ROUND_CEILING: {
449             if (num->sign) needmax=1;   /* Infinity iff non-negative */
450             break;} /* r-c */
451           case DEC_ROUND_FLOOR: {
452             if (!num->sign) needmax=1;  /* Infinity iff negative */
453             break;} /* r-f */
454           default: break;               /* Infinity in all other cases */
455           }
456         if (!needmax) {                 /* easy .. set Infinity */
457           num->exponent=DECFLOAT_Inf;
458           *umsd=0;                      /* be clean: coefficient to 0 */
459           ulsd=umsd;                    /* .. */
460           }
461          else {                         /* return Nmax */
462           umsd=allnines;                /* use constant array */
463           ulsd=allnines+DECPMAX-1;
464           num->exponent=DECEMAX-(DECPMAX-1);
465           }
466         }
467        else { /* no overflow but non-zero and may have to fold-down */
468         Int shift=num->exponent-(DECEMAX-(DECPMAX-1));
469         if (shift>0) {                  /* fold-down needed */
470           /* fold down needed; must copy to buffer in order to pad */
471           /* with zeros safely; fortunately this is not the worst case */
472           /* path because cannot have had a round */
473           uByte buffer[ROUNDUP(DECPMAX+3, 4)]; /* [+3 allows uInt padding] */
474           uByte *s=umsd;                /* source */
475           uByte *t=buffer;              /* safe target */
476           uByte *tlsd=buffer+(ulsd-umsd)+shift; /* target LSD */
477           /* printf("folddown shift=%ld\n", (LI)shift); */
478           for (; s<=ulsd; s+=4, t+=4) UINTAT(t)=UINTAT(s);
479           for (t=tlsd-shift+1; t<=tlsd; t+=4) UINTAT(t)=0;  /* pad */
480           num->exponent-=shift;
481           umsd=buffer;
482           ulsd=tlsd;
483           }
484         } /* fold-down? */
485       length=ulsd-umsd+1;               /* recalculate length */
486       } /* high-end edge case */
487     } /* finite number */
488
489   /*------------------------------------------------------------------*/
490   /* At this point the result will properly fit the decFloat          */
491   /* encoding, and it can be encoded with no possibility of error     */
492   /*------------------------------------------------------------------*/
493   /* Following code does not alter coefficient (could be allnines array) */
494
495   if (length==DECPMAX) {
496     return decFloatFromBCD(df, num->exponent, umsd, num->sign);
497     }
498
499   /* Here when length is short */
500   if (!NUMISSPECIAL(num)) {             /* is still finite */
501     /* encode the combination field and exponent continuation */
502     uInt uexp=(uInt)(num->exponent+DECBIAS); /* biased exponent */
503     uInt code=(uexp>>DECECONL)<<4;      /* top two bits of exp */
504     /* [msd=0] */
505     /* look up the combination field and make high word */
506     encode=DECCOMBFROM[code];           /* indexed by (0-2)*16+msd */
507     encode|=(uexp<<(32-6-DECECONL)) & 0x03ffffff; /* exponent continuation */
508     }
509    else encode=num->exponent;           /* special [already in word] */
510   /* [coefficient length here will be < DECPMAX] */
511
512   encode|=num->sign;                    /* add sign */
513
514   /* private macro to extract a declet, n (where 0<=n<DECLETS and 0 */
515   /* refers to the declet from the least significant three digits) */
516   /* and put the corresponding DPD code into dpd.  Access to umsd and */
517   /* ulsd (pointers to the most and least significant digit of the */
518   /* variable-length coefficient) is assumed, along with use of a */
519   /* working pointer, uInt *ub. */
520   /* As not full-length then chances are there are many leading zeros */
521   /* [and there may be a partial triad] */
522   #define getDPD(dpd, n) ub=ulsd-(3*(n))-2;                           \
523     if (ub<umsd-2) dpd=0;                                             \
524      else if (ub>=umsd) dpd=BCD2DPD[(*ub*256)+(*(ub+1)*16)+*(ub+2)];  \
525      else {dpd=*(ub+2); if (ub+1==umsd) dpd+=*(ub+1)*16; dpd=BCD2DPD[dpd];}
526
527   /* place the declets in the encoding words and copy to result (df), */
528   /* according to endianness; in all cases complete the sign word */
529   /* first */
530   #if DECPMAX==7
531     getDPD(dpd, 1);
532     encode|=dpd<<10;
533     getDPD(dpd, 0);
534     encode|=dpd;
535     DFWORD(df, 0)=encode;     /* just the one word */
536
537   #elif DECPMAX==16
538     getDPD(dpd, 4); encode|=dpd<<8;
539     getDPD(dpd, 3); encode|=dpd>>2;
540     DFWORD(df, 0)=encode;
541     encode=dpd<<30;
542     getDPD(dpd, 2); encode|=dpd<<20;
543     getDPD(dpd, 1); encode|=dpd<<10;
544     getDPD(dpd, 0); encode|=dpd;
545     DFWORD(df, 1)=encode;
546
547   #elif DECPMAX==34
548     getDPD(dpd,10); encode|=dpd<<4;
549     getDPD(dpd, 9); encode|=dpd>>6;
550     DFWORD(df, 0)=encode;
551
552     encode=dpd<<26;
553     getDPD(dpd, 8); encode|=dpd<<16;
554     getDPD(dpd, 7); encode|=dpd<<6;
555     getDPD(dpd, 6); encode|=dpd>>4;
556     DFWORD(df, 1)=encode;
557
558     encode=dpd<<28;
559     getDPD(dpd, 5); encode|=dpd<<18;
560     getDPD(dpd, 4); encode|=dpd<<8;
561     getDPD(dpd, 3); encode|=dpd>>2;
562     DFWORD(df, 2)=encode;
563
564     encode=dpd<<30;
565     getDPD(dpd, 2); encode|=dpd<<20;
566     getDPD(dpd, 1); encode|=dpd<<10;
567     getDPD(dpd, 0); encode|=dpd;
568     DFWORD(df, 3)=encode;
569   #endif
570
571   /* printf("Status: %08lx\n", (LI)set->status); */
572   /* decFloatShow(df, "final"); */
573   return df;
574   } /* decFinalize */
575
576 /* ------------------------------------------------------------------ */
577 /* decFloatFromBCD -- set decFloat from exponent, BCD8, and sign      */
578 /*                                                                    */
579 /*  df is the target decFloat                                         */
580 /*  exp is the in-range unbiased exponent, q, or a special value in   */
581 /*    the form returned by decFloatGetExponent                        */
582 /*  bcdar holds DECPMAX digits to set the coefficient from, one       */
583 /*    digit in each byte (BCD8 encoding); the first (MSD) is ignored  */
584 /*    if df is a NaN; all are ignored if df is infinite.              */
585 /*    All bytes must be in 0-9; results undefined otherwise.          */
586 /*  sig is DECFLOAT_Sign to set the sign bit, 0 otherwise             */
587 /*  returns df, which will be canonical                               */
588 /*                                                                    */
589 /* No error is possible, and no status will be set.                   */
590 /* ------------------------------------------------------------------ */
591 decFloat * decFloatFromBCD(decFloat *df, Int exp, const uByte *bcdar,
592                            Int sig) {
593   uInt encode, dpd;                     /* work */
594   const uByte *ub;                      /* .. */
595
596   if (EXPISSPECIAL(exp)) encode=exp|sig;/* specials already encoded */
597    else {                               /* is finite */
598     /* encode the combination field and exponent continuation */
599     uInt uexp=(uInt)(exp+DECBIAS);      /* biased exponent */
600     uInt code=(uexp>>DECECONL)<<4;      /* top two bits of exp */
601     code+=bcdar[0];                     /* add msd */
602     /* look up the combination field and make high word */
603     encode=DECCOMBFROM[code]|sig;       /* indexed by (0-2)*16+msd */
604     encode|=(uexp<<(32-6-DECECONL)) & 0x03ffffff; /* exponent continuation */
605     }
606
607   /* private macro to extract a declet, n (where 0<=n<DECLETS and 0 */
608   /* refers to the declet from the least significant three digits) */
609   /* and put the corresponding DPD code into dpd. */
610   /* Use of a working pointer, uInt *ub, is assumed. */
611
612   #define getDPDf(dpd, n) ub=bcdar+DECPMAX-1-(3*(n))-2;     \
613     dpd=BCD2DPD[(*ub*256)+(*(ub+1)*16)+*(ub+2)];
614
615   /* place the declets in the encoding words and copy to result (df), */
616   /* according to endianness; in all cases complete the sign word */
617   /* first */
618   #if DECPMAX==7
619     getDPDf(dpd, 1);
620     encode|=dpd<<10;
621     getDPDf(dpd, 0);
622     encode|=dpd;
623     DFWORD(df, 0)=encode;     /* just the one word */
624
625   #elif DECPMAX==16
626     getDPDf(dpd, 4); encode|=dpd<<8;
627     getDPDf(dpd, 3); encode|=dpd>>2;
628     DFWORD(df, 0)=encode;
629     encode=dpd<<30;
630     getDPDf(dpd, 2); encode|=dpd<<20;
631     getDPDf(dpd, 1); encode|=dpd<<10;
632     getDPDf(dpd, 0); encode|=dpd;
633     DFWORD(df, 1)=encode;
634
635   #elif DECPMAX==34
636     getDPDf(dpd,10); encode|=dpd<<4;
637     getDPDf(dpd, 9); encode|=dpd>>6;
638     DFWORD(df, 0)=encode;
639
640     encode=dpd<<26;
641     getDPDf(dpd, 8); encode|=dpd<<16;
642     getDPDf(dpd, 7); encode|=dpd<<6;
643     getDPDf(dpd, 6); encode|=dpd>>4;
644     DFWORD(df, 1)=encode;
645
646     encode=dpd<<28;
647     getDPDf(dpd, 5); encode|=dpd<<18;
648     getDPDf(dpd, 4); encode|=dpd<<8;
649     getDPDf(dpd, 3); encode|=dpd>>2;
650     DFWORD(df, 2)=encode;
651
652     encode=dpd<<30;
653     getDPDf(dpd, 2); encode|=dpd<<20;
654     getDPDf(dpd, 1); encode|=dpd<<10;
655     getDPDf(dpd, 0); encode|=dpd;
656     DFWORD(df, 3)=encode;
657   #endif
658   /* decFloatShow(df, "final"); */
659   return df;
660   } /* decFloatFromBCD */
661
662 /* ------------------------------------------------------------------ */
663 /* decFloatFromPacked -- set decFloat from exponent and packed BCD    */
664 /*                                                                    */
665 /*  df is the target decFloat                                         */
666 /*  exp is the in-range unbiased exponent, q, or a special value in   */
667 /*    the form returned by decFloatGetExponent                        */
668 /*  packed holds DECPMAX packed decimal digits plus a sign nibble     */
669 /*    (all 6 codes are OK); the first (MSD) is ignored if df is a NaN */
670 /*    and all except sign are ignored if df is infinite.  For DOUBLE  */
671 /*    and QUAD the first (pad) nibble is also ignored in all cases.   */
672 /*    All coefficient nibbles must be in 0-9 and sign in A-F; results */
673 /*    are undefined otherwise.                                        */
674 /*  returns df, which will be canonical                               */
675 /*                                                                    */
676 /* No error is possible, and no status will be set.                   */
677 /* ------------------------------------------------------------------ */
678 decFloat * decFloatFromPacked(decFloat *df, Int exp, const uByte *packed) {
679   uByte bcdar[DECPMAX+2];               /* work [+1 for pad, +1 for sign] */
680   const uByte *ip;                      /* .. */
681   uByte *op;                            /* .. */
682   Int   sig=0;                          /* sign */
683
684   /* expand coefficient and sign to BCDAR */
685   #if SINGLE
686   op=bcdar+1;                           /* no pad digit */
687   #else
688   op=bcdar;                             /* first (pad) digit ignored */
689   #endif
690   for (ip=packed; ip<packed+((DECPMAX+2)/2); ip++) {
691     *op++=*ip>>4;
692     *op++=(uByte)(*ip&0x0f);            /* [final nibble is sign] */
693     }
694   op--;                                 /* -> sign byte */
695   if (*op==DECPMINUS || *op==DECPMINUSALT) sig=DECFLOAT_Sign;
696
697   if (EXPISSPECIAL(exp)) {              /* Infinity or NaN */
698     if (!EXPISINF(exp)) bcdar[1]=0;     /* a NaN: ignore MSD */
699      else memset(bcdar+1, 0, DECPMAX);  /* Infinite: coefficient to 0 */
700     }
701   return decFloatFromBCD(df, exp, bcdar+1, sig);
702   } /* decFloatFromPacked */
703
704 /* ------------------------------------------------------------------ */
705 /* decFloatFromString -- conversion from numeric string               */
706 /*                                                                    */
707 /*  result  is the decFloat format number which gets the result of    */
708 /*          the conversion                                            */
709 /*  *string is the character string which should contain a valid      */
710 /*          number (which may be a special value), \0-terminated      */
711 /*          If there are too many significant digits in the           */
712 /*          coefficient it will be rounded.                           */
713 /*  set     is the context                                            */
714 /*  returns result                                                    */
715 /*                                                                    */
716 /* The length of the coefficient and the size of the exponent are     */
717 /* checked by this routine, so the correct error (Underflow or        */
718 /* Overflow) can be reported or rounding applied, as necessary.       */
719 /*                                                                    */
720 /* There is no limit to the coefficient length for finite inputs;     */
721 /* NaN payloads must be integers with no more than DECPMAX-1 digits.  */
722 /* Exponents may have up to nine significant digits.                  */
723 /*                                                                    */
724 /* If bad syntax is detected, the result will be a quiet NaN.         */
725 /* ------------------------------------------------------------------ */
726 decFloat * decFloatFromString(decFloat *result, const char *string,
727                               decContext *set) {
728   Int    digits;                   /* count of digits in coefficient */
729   const  char *dotchar=NULL;       /* where dot was found [NULL if none] */
730   const  char *cfirst=string;      /* -> first character of decimal part */
731   const  char *c;                  /* work */
732   uByte *ub;                       /* .. */
733   bcdnum num;                      /* collects data for finishing */
734   uInt   error=DEC_Conversion_syntax;   /* assume the worst */
735   uByte  buffer[ROUNDUP(DECSTRING+11, 8)]; /* room for most coefficents, */
736                                    /* some common rounding, +3, & pad */
737   #if DECTRACE
738   /* printf("FromString %s ...\n", string); */
739   #endif
740
741   for(;;) {                             /* once-only 'loop' */
742     num.sign=0;                         /* assume non-negative */
743     num.msd=buffer;                     /* MSD is here always */
744
745     /* detect and validate the coefficient, including any leading, */
746     /* trailing, or embedded '.' */
747     /* [could test four-at-a-time here (saving 10% for decQuads), */
748     /* but that risks storage violation because the position of the */
749     /* terminator is unknown] */
750     for (c=string;; c++) {              /* -> input character */
751       if (((unsigned)(*c-'0'))<=9) continue; /* '0' through '9' is good */
752       if (*c=='\0') break;              /* most common non-digit */
753       if (*c=='.') {
754         if (dotchar!=NULL) break;       /* not first '.' */
755         dotchar=c;                      /* record offset into decimal part */
756         continue;}
757       if (c==string) {                  /* first in string... */
758         if (*c=='-') {                  /* valid - sign */
759           cfirst++;
760           num.sign=DECFLOAT_Sign;
761           continue;}
762         if (*c=='+') {                  /* valid + sign */
763           cfirst++;
764           continue;}
765         }
766       /* *c is not a digit, terminator, or a valid +, -, or '.' */
767       break;
768       } /* c loop */
769
770     digits=(uInt)(c-cfirst);            /* digits (+1 if a dot) */
771
772     if (digits>0) {                     /* had digits and/or dot */
773       const char *clast=c-1;            /* note last coefficient char position */
774       Int exp=0;                        /* exponent accumulator */
775       if (*c!='\0') {                   /* something follows the coefficient */
776         uInt edig;                      /* unsigned work */
777         /* had some digits and more to come; expect E[+|-]nnn now */
778         const char *firstexp;           /* exponent first non-zero */
779         if (*c!='E' && *c!='e') break;
780         c++;                            /* to (optional) sign */
781         if (*c=='-' || *c=='+') c++;    /* step over sign (c=clast+2) */
782         if (*c=='\0') break;            /* no digits!  (e.g., '1.2E') */
783         for (; *c=='0';) c++;           /* skip leading zeros [even last] */
784         firstexp=c;                     /* remember start [maybe '\0'] */
785         /* gather exponent digits */
786         edig=(uInt)*c-(uInt)'0';
787         if (edig<=9) {                  /* [check not bad or terminator] */
788           exp+=edig;                    /* avoid initial X10 */
789           c++;
790           for (;; c++) {
791             edig=(uInt)*c-(uInt)'0';
792             if (edig>9) break;
793             exp=exp*10+edig;
794             }
795           }
796         /* if not now on the '\0', *c must not be a digit */
797         if (*c!='\0') break;
798
799         /* (this next test must be after the syntax checks) */
800         /* if definitely more than the possible digits for format then */
801         /* the exponent may have wrapped, so simply set it to a certain */
802         /* over/underflow value */
803         if (c>firstexp+DECEMAXD) exp=DECEMAX*2;
804         if (*(clast+2)=='-') exp=-exp;  /* was negative */
805         } /* digits>0 */
806
807       if (dotchar!=NULL) {              /* had a '.' */
808         digits--;                       /* remove from digits count */
809         if (digits==0) break;           /* was dot alone: bad syntax */
810         exp-=(Int)(clast-dotchar);      /* adjust exponent */
811         /* [the '.' can now be ignored] */
812         }
813       num.exponent=exp;                 /* exponent is good; store it */
814
815       /* Here when whole string has been inspected and syntax is good */
816       /* cfirst->first digit or dot, clast->last digit or dot */
817       error=0;                          /* no error possible now */
818
819       /* if the number of digits in the coefficient will fit in buffer */
820       /* then it can simply be converted to bcd8 and copied -- decFinalize */
821       /* will take care of leading zeros and rounding; the buffer is big */
822       /* enough for all canonical coefficients, including 0.00000nn... */
823       ub=buffer;
824       if (digits<=(Int)(sizeof(buffer)-3)) { /* [-3 allows by-4s copy] */
825         c=cfirst;
826         if (dotchar!=NULL) {                 /* a dot to worry about */
827           if (*(c+1)=='.') {                 /* common canonical case */
828             *ub++=(uByte)(*c-'0');           /* copy leading digit */
829             c+=2;                            /* prepare to handle rest */
830             }
831            else for (; c<=clast;) {          /* '.' could be anywhere */
832             /* as usual, go by fours when safe; NB it has been asserted */
833             /* that a '.' does not have the same mask as a digit */
834             if (c<=clast-3                             /* safe for four */
835              && (UINTAT(c)&0xf0f0f0f0)==CHARMASK) {    /* test four */
836               UINTAT(ub)=UINTAT(c)&0x0f0f0f0f;         /* to BCD8 */
837               ub+=4;
838               c+=4;
839               continue;
840               }
841             if (*c=='.') {                   /* found the dot */
842               c++;                           /* step over it .. */
843               break;                         /* .. and handle the rest */
844               }
845             *ub++=(uByte)(*c++-'0');
846             }
847           } /* had dot */
848         /* Now no dot; do this by fours (where safe) */
849         for (; c<=clast-3; c+=4, ub+=4) UINTAT(ub)=UINTAT(c)&0x0f0f0f0f;
850         for (; c<=clast; c++, ub++) *ub=(uByte)(*c-'0');
851         num.lsd=buffer+digits-1;             /* record new LSD */
852         } /* fits */
853
854        else {                                /* too long for buffer */
855         /* [This is a rare and unusual case; arbitrary-length input] */
856         /* strip leading zeros [but leave final 0 if all 0's] */
857         if (*cfirst=='.') cfirst++;          /* step past dot at start */
858         if (*cfirst=='0') {                  /* [cfirst always -> digit] */
859           for (; cfirst<clast; cfirst++) {
860             if (*cfirst!='0') {              /* non-zero found */
861               if (*cfirst=='.') continue;    /* [ignore] */
862               break;                         /* done */
863               }
864             digits--;                        /* 0 stripped */
865             } /* cfirst */
866           } /* at least one leading 0 */
867
868         /* the coefficient is now as short as possible, but may still */
869         /* be too long; copy up to Pmax+1 digits to the buffer, then */
870         /* just record any non-zeros (set round-for-reround digit) */
871         for (c=cfirst; c<=clast && ub<=buffer+DECPMAX; c++) {
872           /* (see commentary just above) */
873           if (c<=clast-3                          /* safe for four */
874            && (UINTAT(c)&0xf0f0f0f0)==CHARMASK) { /* four digits */
875             UINTAT(ub)=UINTAT(c)&0x0f0f0f0f;      /* to BCD8 */
876             ub+=4;
877             c+=3;                            /* [will become 4] */
878             continue;
879             }
880           if (*c=='.') continue;             /* [ignore] */
881           *ub++=(uByte)(*c-'0');
882           }
883         ub--;                                /* -> LSD */
884         for (; c<=clast; c++) {              /* inspect remaining chars */
885           if (*c!='0') {                     /* sticky bit needed */
886             if (*c=='.') continue;           /* [ignore] */
887             *ub=DECSTICKYTAB[*ub];           /* update round-for-reround */
888             break;                           /* no need to look at more */
889             }
890           }
891         num.lsd=ub;                          /* record LSD */
892         /* adjust exponent for dropped digits */
893         num.exponent+=digits-(Int)(ub-buffer+1);
894         } /* too long for buffer */
895       } /* digits or dot */
896
897      else {                             /* no digits or dot were found */
898       if (*c=='\0') break;              /* nothing to come is bad */
899       /* only Infinities and NaNs are allowed, here */
900       buffer[0]=0;                      /* default a coefficient of 0 */
901       num.lsd=buffer;                   /* .. */
902       if (decBiStr(c, "infinity", "INFINITY")
903        || decBiStr(c, "inf", "INF")) num.exponent=DECFLOAT_Inf;
904        else {                           /* should be a NaN */
905         num.exponent=DECFLOAT_qNaN;     /* assume quiet NaN */
906         if (*c=='s' || *c=='S') {       /* probably an sNaN */
907           c++;
908           num.exponent=DECFLOAT_sNaN;   /* assume is in fact sNaN */
909           }
910         if (*c!='N' && *c!='n') break;  /* check caseless "NaN" */
911         c++;
912         if (*c!='a' && *c!='A') break;  /* .. */
913         c++;
914         if (*c!='N' && *c!='n') break;  /* .. */
915         c++;
916         /* now either nothing, or nnnn payload (no dots), expected */
917         /* -> start of integer, and skip leading 0s [including plain 0] */
918         for (cfirst=c; *cfirst=='0';) cfirst++;
919         if (*cfirst!='\0') {            /* not empty or all-0, payload */
920           /* payload found; check all valid digits and copy to buffer as bcd8 */
921           ub=buffer;
922           for (c=cfirst;; c++, ub++) {
923             if ((unsigned)(*c-'0')>9) break; /* quit if not 0-9 */
924             if (c-cfirst==DECPMAX-1) break;  /* too many digits */
925             *ub=(uByte)(*c-'0');        /* good bcd8 */
926             }
927           if (*c!='\0') break;          /* not all digits, or too many */
928           num.lsd=ub-1;                 /* record new LSD */
929           }
930         } /* NaN or sNaN */
931       error=0;                          /* syntax is OK */
932       break;                            /* done with specials */
933       } /* digits=0 (special expected) */
934     break;
935     }                                   /* [for(;;) break] */
936
937   /* decShowNum(&num, "fromStr"); */
938
939   if (error!=0) {
940     set->status|=error;
941     num.exponent=DECFLOAT_qNaN;         /* set up quiet NaN */
942     num.sign=0;                         /* .. with 0 sign */
943     buffer[0]=0;                        /* .. and coefficient */
944     num.lsd=buffer;                     /* .. */
945     /* decShowNum(&num, "oops"); */
946     }
947
948   /* decShowNum(&num, "dffs"); */
949   decFinalize(result, &num, set);       /* round, check, and lay out */
950   /* decFloatShow(result, "fromString"); */
951   return result;
952   } /* decFloatFromString */
953
954 /* ------------------------------------------------------------------ */
955 /* decFloatFromWider -- conversion from next-wider format             */
956 /*                                                                    */
957 /*  result  is the decFloat format number which gets the result of    */
958 /*          the conversion                                            */
959 /*  wider   is the decFloatWider format number which will be narrowed */
960 /*  set     is the context                                            */
961 /*  returns result                                                    */
962 /*                                                                    */
963 /* Narrowing can cause rounding, overflow, etc., but not Invalid      */
964 /* operation (sNaNs are copied and do not signal).                    */
965 /* ------------------------------------------------------------------ */
966 /* narrow-to is not possible for decQuad format numbers; simply omit */
967 #if !QUAD
968 decFloat * decFloatFromWider(decFloat *result, const decFloatWider *wider,
969                              decContext *set) {
970   bcdnum num;                           /* collects data for finishing */
971   uByte  bcdar[DECWPMAX];               /* room for wider coefficient */
972   uInt   widerhi=DFWWORD(wider, 0);     /* top word */
973   Int    exp;
974
975   GETWCOEFF(wider, bcdar);
976
977   num.msd=bcdar;                        /* MSD is here always */
978   num.lsd=bcdar+DECWPMAX-1;             /* LSD is here always */
979   num.sign=widerhi&0x80000000;          /* extract sign [DECFLOAT_Sign=Neg] */
980
981   /* decode the wider combination field to exponent */
982   exp=DECCOMBWEXP[widerhi>>26];         /* decode from wider combination field */
983   /* if it is a special there's nothing to do unless sNaN; if it's */
984   /* finite then add the (wider) exponent continuation and unbias */
985   if (EXPISSPECIAL(exp)) exp=widerhi&0x7e000000; /* include sNaN selector */
986    else exp+=GETWECON(wider)-DECWBIAS;
987   num.exponent=exp;
988
989   /* decShowNum(&num, "dffw"); */
990   return decFinalize(result, &num, set);/* round, check, and lay out */
991   } /* decFloatFromWider */
992 #endif
993
994 /* ------------------------------------------------------------------ */
995 /* decFloatGetCoefficient -- get coefficient as BCD8                  */
996 /*                                                                    */
997 /*  df is the decFloat from which to extract the coefficient          */
998 /*  bcdar is where DECPMAX bytes will be written, one BCD digit in    */
999 /*    each byte (BCD8 encoding); if df is a NaN the first byte will   */
1000 /*    be zero, and if it is infinite they will all be zero            */
1001 /*  returns the sign of the coefficient (DECFLOAT_Sign if negative,   */
1002 /*    0 otherwise)                                                    */
1003 /*                                                                    */
1004 /* No error is possible, and no status will be set.  If df is a       */
1005 /* special value the array is set to zeros (for Infinity) or to the   */
1006 /* payload of a qNaN or sNaN.                                         */
1007 /* ------------------------------------------------------------------ */
1008 Int decFloatGetCoefficient(const decFloat *df, uByte *bcdar) {
1009   if (DFISINF(df)) memset(bcdar, 0, DECPMAX);
1010    else {
1011     GETCOEFF(df, bcdar);           /* use macro */
1012     if (DFISNAN(df)) bcdar[0]=0;   /* MSD needs correcting */
1013     }
1014   return DFISSIGNED(df);
1015   } /* decFloatGetCoefficient */
1016
1017 /* ------------------------------------------------------------------ */
1018 /* decFloatGetExponent -- get unbiased exponent                       */
1019 /*                                                                    */
1020 /*  df is the decFloat from which to extract the exponent             */
1021 /*  returns the exponent, q.                                          */
1022 /*                                                                    */
1023 /* No error is possible, and no status will be set.  If df is a       */
1024 /* special value the first seven bits of the decFloat are returned,   */
1025 /* left adjusted and with the first (sign) bit set to 0 (followed by  */
1026 /* 25 0 bits).  e.g., -sNaN would return 0x7e000000 (DECFLOAT_sNaN).  */
1027 /* ------------------------------------------------------------------ */
1028 Int decFloatGetExponent(const decFloat *df) {
1029   if (DFISSPECIAL(df)) return DFWORD(df, 0)&0x7e000000;
1030   return GETEXPUN(df);
1031   } /* decFloatGetExponent */
1032
1033 /* ------------------------------------------------------------------ */
1034 /* decFloatSetCoefficient -- set coefficient from BCD8                */
1035 /*                                                                    */
1036 /*  df is the target decFloat (and source of exponent/special value)  */
1037 /*  bcdar holds DECPMAX digits to set the coefficient from, one       */
1038 /*    digit in each byte (BCD8 encoding); the first (MSD) is ignored  */
1039 /*    if df is a NaN; all are ignored if df is infinite.              */
1040 /*  sig is DECFLOAT_Sign to set the sign bit, 0 otherwise             */
1041 /*  returns df, which will be canonical                               */
1042 /*                                                                    */
1043 /* No error is possible, and no status will be set.                   */
1044 /* ------------------------------------------------------------------ */
1045 decFloat * decFloatSetCoefficient(decFloat *df, const uByte *bcdar,
1046                                   Int sig) {
1047   uInt exp;                        /* for exponent */
1048   uByte bcdzero[DECPMAX];          /* for infinities */
1049
1050   /* Exponent/special code is extracted from df */
1051   if (DFISSPECIAL(df)) {
1052     exp=DFWORD(df, 0)&0x7e000000;
1053     if (DFISINF(df)) {
1054       memset(bcdzero, 0, DECPMAX);
1055       return decFloatFromBCD(df, exp, bcdzero, sig);
1056       }
1057     }
1058    else exp=GETEXPUN(df);
1059   return decFloatFromBCD(df, exp, bcdar, sig);
1060   } /* decFloatSetCoefficient */
1061
1062 /* ------------------------------------------------------------------ */
1063 /* decFloatSetExponent -- set exponent or special value               */
1064 /*                                                                    */
1065 /*  df  is the target decFloat (and source of coefficient/payload)    */
1066 /*  set is the context for reporting status                           */
1067 /*  exp is the unbiased exponent, q, or a special value in the form   */
1068 /*    returned by decFloatGetExponent                                 */
1069 /*  returns df, which will be canonical                               */
1070 /*                                                                    */
1071 /* No error is possible, but Overflow or Underflow might occur.       */
1072 /* ------------------------------------------------------------------ */
1073 decFloat * decFloatSetExponent(decFloat *df, decContext *set, Int exp) {
1074   uByte  bcdcopy[DECPMAX];         /* for coefficient */
1075   bcdnum num;                      /* work */
1076   num.exponent=exp;
1077   num.sign=decFloatGetCoefficient(df, bcdcopy); /* extract coefficient */
1078   if (DFISSPECIAL(df)) {           /* MSD or more needs correcting */
1079     if (DFISINF(df)) memset(bcdcopy, 0, DECPMAX);
1080     bcdcopy[0]=0;
1081     }
1082   num.msd=bcdcopy;
1083   num.lsd=bcdcopy+DECPMAX-1;
1084   return decFinalize(df, &num, set);
1085   } /* decFloatSetExponent */
1086
1087 /* ------------------------------------------------------------------ */
1088 /* decFloatRadix -- returns the base (10)                             */
1089 /*                                                                    */
1090 /*   df is any decFloat of this format                                */
1091 /* ------------------------------------------------------------------ */
1092 uInt decFloatRadix(const decFloat *df) {
1093   if (df) return 10;                         /* to placate compiler */
1094   return 10;
1095   } /* decFloatRadix */
1096
1097 /* ------------------------------------------------------------------ */
1098 /* decFloatShow -- printf a decFloat in hexadecimal and decimal       */
1099 /*   df  is the decFloat to show                                      */
1100 /*   tag is a tag string displayed with the number                    */
1101 /*                                                                    */
1102 /* This is a debug aid; the precise format of the string may change.  */
1103 /* ------------------------------------------------------------------ */
1104 void decFloatShow(const decFloat *df, const char *tag) {
1105   char hexbuf[DECBYTES*2+DECBYTES/4+1]; /* NB blank after every fourth */
1106   char buff[DECSTRING];                 /* for value in decimal */
1107   Int i, j=0;
1108
1109   for (i=0; i<DECBYTES; i++) {
1110     #if DECLITEND
1111       sprintf(&hexbuf[j], "%02x", df->bytes[DECBYTES-1-i]);
1112     #else
1113       sprintf(&hexbuf[j], "%02x", df->bytes[i]);
1114     #endif
1115     j+=2;
1116     /* the next line adds blank (and terminator) after final pair, too */
1117     if ((i+1)%4==0) {strcpy(&hexbuf[j], " "); j++;}
1118     }
1119   decFloatToString(df, buff);
1120   printf(">%s> %s [big-endian]  %s\n", tag, hexbuf, buff);
1121   return;
1122   } /* decFloatShow */
1123
1124 /* ------------------------------------------------------------------ */
1125 /* decFloatToBCD -- get sign, exponent, and BCD8 from a decFloat      */
1126 /*                                                                    */
1127 /*  df is the source decFloat                                         */
1128 /*  exp will be set to the unbiased exponent, q, or to a special      */
1129 /*    value in the form returned by decFloatGetExponent               */
1130 /*  bcdar is where DECPMAX bytes will be written, one BCD digit in    */
1131 /*    each byte (BCD8 encoding); if df is a NaN the first byte will   */
1132 /*    be zero, and if it is infinite they will all be zero            */
1133 /*  returns the sign of the coefficient (DECFLOAT_Sign if negative,   */
1134 /*    0 otherwise)                                                    */
1135 /*                                                                    */
1136 /* No error is possible, and no status will be set.                   */
1137 /* ------------------------------------------------------------------ */
1138 Int decFloatToBCD(const decFloat *df, Int *exp, uByte *bcdar) {
1139   if (DFISINF(df)) {
1140     memset(bcdar, 0, DECPMAX);
1141     *exp=DFWORD(df, 0)&0x7e000000;
1142     }
1143    else {
1144     GETCOEFF(df, bcdar);           /* use macro */
1145     if (DFISNAN(df)) {
1146       bcdar[0]=0;                  /* MSD needs correcting */
1147       *exp=DFWORD(df, 0)&0x7e000000;
1148       }
1149      else {                        /* finite */
1150       *exp=GETEXPUN(df);
1151       }
1152     }
1153   return DFISSIGNED(df);
1154   } /* decFloatToBCD */
1155
1156 /* ------------------------------------------------------------------ */
1157 /* decFloatToEngString -- conversion to numeric string, engineering   */
1158 /*                                                                    */
1159 /*  df is the decFloat format number to convert                       */
1160 /*  string is the string where the result will be laid out            */
1161 /*                                                                    */
1162 /* string must be at least DECPMAX+9 characters (the worst case is    */
1163 /* "-0.00000nnn...nnn\0", which is as long as the exponent form when  */
1164 /* DECEMAXD<=4); this condition is asserted above                     */
1165 /*                                                                    */
1166 /* No error is possible, and no status will be set                    */
1167 /* ------------------------------------------------------------------ */
1168 char * decFloatToEngString(const decFloat *df, char *string){
1169   uInt msd;                        /* coefficient MSD */
1170   Int  exp;                        /* exponent top two bits or full */
1171   uInt comb;                       /* combination field */
1172   char *cstart;                    /* coefficient start */
1173   char *c;                         /* output pointer in string */
1174   char *s, *t;                     /* .. (source, target) */
1175   Int  pre, e;                     /* work */
1176   const uByte *u;                  /* .. */
1177
1178   /* Source words; macro handles endianness */
1179   uInt sourhi=DFWORD(df, 0);       /* word with sign */
1180   #if DECPMAX==16
1181   uInt sourlo=DFWORD(df, 1);
1182   #elif DECPMAX==34
1183   uInt sourmh=DFWORD(df, 1);
1184   uInt sourml=DFWORD(df, 2);
1185   uInt sourlo=DFWORD(df, 3);
1186   #endif
1187
1188   c=string;                        /* where result will go */
1189   if (((Int)sourhi)<0) *c++='-';   /* handle sign */
1190   comb=sourhi>>26;                 /* sign+combination field */
1191   msd=DECCOMBMSD[comb];            /* decode the combination field */
1192   exp=DECCOMBEXP[comb];            /* .. */
1193
1194   if (EXPISSPECIAL(exp)) {         /* special */
1195     if (exp==DECFLOAT_Inf) {       /* infinity */
1196       strcpy(c,   "Inf");
1197       strcpy(c+3, "inity");
1198       return string;               /* easy */
1199       }
1200     if (sourhi&0x02000000) *c++='s'; /* sNaN */
1201     strcpy(c, "NaN");              /* complete word */
1202     c+=3;                          /* step past */
1203     /* quick exit if the payload is zero */
1204     #if DECPMAX==7
1205     if ((sourhi&0x000fffff)==0) return string;
1206     #elif DECPMAX==16
1207     if (sourlo==0 && (sourhi&0x0003ffff)==0) return string;
1208     #elif DECPMAX==34
1209     if (sourlo==0 && sourml==0 && sourmh==0
1210      && (sourhi&0x00003fff)==0) return string;
1211     #endif
1212     /* otherwise drop through to add integer; set correct exp etc. */
1213     exp=0; msd=0;                  /* setup for following code */
1214     }
1215    else { /* complete exponent; top two bits are in place */
1216     exp+=GETECON(df)-DECBIAS;      /* .. + continuation and unbias */
1217     }
1218
1219   /* convert the digits of the significand to characters */
1220   cstart=c;                        /* save start of coefficient */
1221   if (msd) *c++=(char)('0'+(char)msd);  /* non-zero most significant digit */
1222
1223   /* Decode the declets.  After extracting each declet, it is */
1224   /* decoded to a 4-uByte sequence by table lookup; the four uBytes */
1225   /* are the three encoded BCD8 digits followed by a 1-byte length */
1226   /* (significant digits, except that 000 has length 0).  This allows */
1227   /* us to left-align the first declet with non-zero content, then */
1228   /* the remaining ones are full 3-char length.  Fixed-length copies */
1229   /* are used because variable-length memcpy causes a subroutine call */
1230   /* in at least two compilers.  (The copies are length 4 for speed */
1231   /* and are safe because the last item in the array is of length */
1232   /* three and has the length byte following.) */
1233   #define dpd2char(dpdin) u=&DPD2BCD8[((dpdin)&0x3ff)*4];        \
1234          if (c!=cstart) {UINTAT(c)=UINTAT(u)|CHARMASK; c+=3;}    \
1235           else if (*(u+3)) {                                     \
1236            UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; c+=*(u+3);}
1237
1238   #if DECPMAX==7
1239   dpd2char(sourhi>>10);                 /* declet 1 */
1240   dpd2char(sourhi);                     /* declet 2 */
1241
1242   #elif DECPMAX==16
1243   dpd2char(sourhi>>8);                  /* declet 1 */
1244   dpd2char((sourhi<<2) | (sourlo>>30)); /* declet 2 */
1245   dpd2char(sourlo>>20);                 /* declet 3 */
1246   dpd2char(sourlo>>10);                 /* declet 4 */
1247   dpd2char(sourlo);                     /* declet 5 */
1248
1249   #elif DECPMAX==34
1250   dpd2char(sourhi>>4);                  /* declet 1 */
1251   dpd2char((sourhi<<6) | (sourmh>>26)); /* declet 2 */
1252   dpd2char(sourmh>>16);                 /* declet 3 */
1253   dpd2char(sourmh>>6);                  /* declet 4 */
1254   dpd2char((sourmh<<4) | (sourml>>28)); /* declet 5 */
1255   dpd2char(sourml>>18);                 /* declet 6 */
1256   dpd2char(sourml>>8);                  /* declet 7 */
1257   dpd2char((sourml<<2) | (sourlo>>30)); /* declet 8 */
1258   dpd2char(sourlo>>20);                 /* declet 9 */
1259   dpd2char(sourlo>>10);                 /* declet 10 */
1260   dpd2char(sourlo);                     /* declet 11 */
1261   #endif
1262
1263   if (c==cstart) *c++='0';         /* all zeros, empty -- make "0" */
1264
1265   if (exp==0) {                    /* integer or NaN case -- easy */
1266     *c='\0';                       /* terminate */
1267     return string;
1268     }
1269   /* non-0 exponent */
1270
1271   e=0;                             /* assume no E */
1272   pre=(Int)(c-cstart)+exp;         /* length+exp  [c->LSD+1] */
1273   /* [here, pre-exp is the digits count (==1 for zero)] */
1274
1275   if (exp>0 || pre<-5) {           /* need exponential form */
1276     e=pre-1;                       /* calculate E value */
1277     pre=1;                         /* assume one digit before '.' */
1278     if (e!=0) {                    /* engineering: may need to adjust */
1279       Int adj;                     /* adjustment */
1280       /* The C remainder operator is undefined for negative numbers, so */
1281       /* a positive remainder calculation must be used here */
1282       if (e<0) {
1283         adj=(-e)%3;
1284         if (adj!=0) adj=3-adj;
1285         }
1286        else { /* e>0 */
1287         adj=e%3;
1288         }
1289       e=e-adj;
1290       /* if dealing with zero still produce an exponent which is a */
1291       /* multiple of three, as expected, but there will only be the */
1292       /* one zero before the E, still.  Otherwise note the padding. */
1293       if (!DFISZERO(df)) pre+=adj;
1294        else {  /* is zero */
1295         if (adj!=0) {              /* 0.00Esnn needed */
1296           e=e+3;
1297           pre=-(2-adj);
1298           }
1299         } /* zero */
1300       } /* engineering adjustment */
1301     } /* exponential form */
1302   /* printf("e=%ld pre=%ld exp=%ld\n", (LI)e, (LI)pre, (LI)exp); */
1303
1304   /* modify the coefficient, adding 0s, '.', and E+nn as needed */
1305   if (pre>0) {                     /* ddd.ddd (plain), perhaps with E */
1306                                    /* or dd00 padding for engineering */
1307     char *dotat=cstart+pre;
1308     if (dotat<c) {                      /* if embedded dot needed... */
1309       /* move by fours; there must be space for junk at the end */
1310       /* because there is still space for exponent */
1311       s=dotat+ROUNDDOWN4(c-dotat);      /* source */
1312       t=s+1;                            /* target */
1313       /* open the gap */
1314       for (; s>=dotat; s-=4, t-=4) UINTAT(t)=UINTAT(s);
1315       *dotat='.';
1316       c++;                              /* length increased by one */
1317       } /* need dot? */
1318      else for (; c<dotat; c++) *c='0';  /* pad for engineering */
1319     } /* pre>0 */
1320    else {
1321     /* -5<=pre<=0: here for plain 0.ddd or 0.000ddd forms (may have
1322        E, but only for 0.00E+3 kind of case -- with plenty of spare
1323        space in this case */
1324     pre=-pre+2;                         /* gap width, including "0." */
1325     t=cstart+ROUNDDOWN4(c-cstart)+pre;  /* preferred first target point */
1326     /* backoff if too far to the right */
1327     if (t>string+DECSTRING-5) t=string+DECSTRING-5; /* adjust to fit */
1328     /* now shift the entire coefficient to the right, being careful not */
1329     /* to access to the left of string */
1330     for (s=t-pre; s>=string; s-=4, t-=4) UINTAT(t)=UINTAT(s);
1331     /* for Quads and Singles there may be a character or two left... */
1332     s+=3;                               /* where next would come from */
1333     for(; s>=cstart; s--, t--) *(t+3)=*(s);
1334     /* now have fill 0. through 0.00000; use overlaps to avoid tests */
1335     if (pre>=4) {
1336       UINTAT(cstart+pre-4)=UINTAT("0000");
1337       UINTAT(cstart)=UINTAT("0.00");
1338       }
1339      else { /* 2 or 3 */
1340       *(cstart+pre-1)='0';
1341       USHORTAT(cstart)=USHORTAT("0.");
1342       }
1343     c+=pre;                             /* to end */
1344     }
1345
1346   /* finally add the E-part, if needed; it will never be 0, and has */
1347   /* a maximum length of 3 or 4 digits (asserted above) */
1348   if (e!=0) {
1349     USHORTAT(c)=USHORTAT("E+");         /* starts with E, assume + */
1350     c++;
1351     if (e<0) {
1352       *c='-';                           /* oops, need '-' */
1353       e=-e;                             /* uInt, please */
1354       }
1355     c++;
1356     /* Three-character exponents are easy; 4-character a little trickier */
1357     #if DECEMAXD<=3
1358       u=&BIN2BCD8[e*4];                 /* -> 3 digits + length byte */
1359       /* copy fixed 4 characters [is safe], starting at non-zero */
1360       /* and with character mask to convert BCD to char */
1361       UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK;
1362       c+=*(u+3);                        /* bump pointer appropriately */
1363     #elif DECEMAXD==4
1364       if (e<1000) {                     /* 3 (or fewer) digits case */
1365         u=&BIN2BCD8[e*4];               /* -> 3 digits + length byte */
1366         UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; /* [as above] */
1367         c+=*(u+3);                      /* bump pointer appropriately */
1368         }
1369        else {                           /* 4-digits */
1370         Int thou=((e>>3)*1049)>>17;     /* e/1000 */
1371         Int rem=e-(1000*thou);          /* e%1000 */
1372         *c++=(char)('0'+(char)thou);    /* the thousands digit */
1373         u=&BIN2BCD8[rem*4];             /* -> 3 digits + length byte */
1374         UINTAT(c)=UINTAT(u)|CHARMASK;   /* copy fixed 3+1 characters [is safe] */
1375         c+=3;                           /* bump pointer, always 3 digits */
1376         }
1377     #endif
1378     }
1379   *c='\0';                              /* terminate */
1380   /*printf("res %s\n", string); */
1381   return string;
1382   } /* decFloatToEngString */
1383
1384 /* ------------------------------------------------------------------ */
1385 /* decFloatToPacked -- convert decFloat to Packed decimal + exponent  */
1386 /*                                                                    */
1387 /*  df is the source decFloat                                         */
1388 /*  exp will be set to the unbiased exponent, q, or to a special      */
1389 /*    value in the form returned by decFloatGetExponent               */
1390 /*  packed is where DECPMAX nibbles will be written with the sign as  */
1391 /*    final nibble (0x0c for +, 0x0d for -); a NaN has a first nibble */
1392 /*    of zero, and an infinity is all zeros. decDouble and decQuad    */
1393 /*    have a additional leading zero nibble, leading to result        */
1394 /*    lengths of 4, 9, and 18 bytes.                                  */
1395 /*  returns the sign of the coefficient (DECFLOAT_Sign if negative,   */
1396 /*    0 otherwise)                                                    */
1397 /*                                                                    */
1398 /* No error is possible, and no status will be set.                   */
1399 /* ------------------------------------------------------------------ */
1400 Int decFloatToPacked(const decFloat *df, Int *exp, uByte *packed) {
1401   uByte bcdar[DECPMAX+2];          /* work buffer */
1402   uByte *ip=bcdar, *op=packed;     /* work pointers */
1403   if (DFISINF(df)) {
1404     memset(bcdar, 0, DECPMAX+2);
1405     *exp=DECFLOAT_Inf;
1406     }
1407    else {
1408     GETCOEFF(df, bcdar+1);         /* use macro */
1409     if (DFISNAN(df)) {
1410       bcdar[1]=0;                  /* MSD needs clearing */
1411       *exp=DFWORD(df, 0)&0x7e000000;
1412       }
1413      else {                        /* finite */
1414       *exp=GETEXPUN(df);
1415       }
1416     }
1417   /* now pack; coefficient currently at bcdar+1 */
1418   #if SINGLE
1419     ip++;                          /* ignore first byte */
1420   #else
1421     *ip=0;                         /* need leading zero */
1422   #endif
1423   /* set final byte to Packed BCD sign value */
1424   bcdar[DECPMAX+1]=(DFISSIGNED(df) ? DECPMINUS : DECPPLUS);
1425   /* pack an even number of bytes... */
1426   for (; op<packed+((DECPMAX+2)/2); op++, ip+=2) {
1427     *op=(uByte)((*ip<<4)+*(ip+1));
1428     }
1429   return (bcdar[DECPMAX+1]==DECPMINUS ? DECFLOAT_Sign : 0);
1430   } /* decFloatToPacked */
1431
1432 /* ------------------------------------------------------------------ */
1433 /* decFloatToString -- conversion to numeric string                   */
1434 /*                                                                    */
1435 /*  df is the decFloat format number to convert                       */
1436 /*  string is the string where the result will be laid out            */
1437 /*                                                                    */
1438 /* string must be at least DECPMAX+9 characters (the worst case is    */
1439 /* "-0.00000nnn...nnn\0", which is as long as the exponent form when  */
1440 /* DECEMAXD<=4); this condition is asserted above                     */
1441 /*                                                                    */
1442 /* No error is possible, and no status will be set                    */
1443 /* ------------------------------------------------------------------ */
1444 char * decFloatToString(const decFloat *df, char *string){
1445   uInt msd;                        /* coefficient MSD */
1446   Int  exp;                        /* exponent top two bits or full */
1447   uInt comb;                       /* combination field */
1448   char *cstart;                    /* coefficient start */
1449   char *c;                         /* output pointer in string */
1450   char *s, *t;                     /* .. (source, target) */
1451   Int  pre, e;                     /* work */
1452   const uByte *u;                  /* .. */
1453
1454   /* Source words; macro handles endianness */
1455   uInt sourhi=DFWORD(df, 0);       /* word with sign */
1456   #if DECPMAX==16
1457   uInt sourlo=DFWORD(df, 1);
1458   #elif DECPMAX==34
1459   uInt sourmh=DFWORD(df, 1);
1460   uInt sourml=DFWORD(df, 2);
1461   uInt sourlo=DFWORD(df, 3);
1462   #endif
1463
1464   c=string;                        /* where result will go */
1465   if (((Int)sourhi)<0) *c++='-';   /* handle sign */
1466   comb=sourhi>>26;                 /* sign+combination field */
1467   msd=DECCOMBMSD[comb];            /* decode the combination field */
1468   exp=DECCOMBEXP[comb];            /* .. */
1469
1470   if (EXPISSPECIAL(exp)) {         /* special */
1471     if (exp==DECFLOAT_Inf) {       /* infinity */
1472       strcpy(c, "Infinity");
1473       return string;               /* easy */
1474       }
1475     if (sourhi&0x02000000) *c++='s'; /* sNaN */
1476     strcpy(c, "NaN");              /* complete word */
1477     c+=3;                          /* step past */
1478     /* quick exit if the payload is zero */
1479     #if DECPMAX==7
1480     if ((sourhi&0x000fffff)==0) return string;
1481     #elif DECPMAX==16
1482     if (sourlo==0 && (sourhi&0x0003ffff)==0) return string;
1483     #elif DECPMAX==34
1484     if (sourlo==0 && sourml==0 && sourmh==0
1485      && (sourhi&0x00003fff)==0) return string;
1486     #endif
1487     /* otherwise drop through to add integer; set correct exp etc. */
1488     exp=0; msd=0;                  /* setup for following code */
1489     }
1490    else { /* complete exponent; top two bits are in place */
1491     exp+=GETECON(df)-DECBIAS;      /* .. + continuation and unbias */
1492     }
1493
1494   /* convert the digits of the significand to characters */
1495   cstart=c;                        /* save start of coefficient */
1496   if (msd) *c++=(char)('0'+(char)msd);  /* non-zero most significant digit */
1497
1498   /* Decode the declets.  After extracting each declet, it is */
1499   /* decoded to a 4-uByte sequence by table lookup; the four uBytes */
1500   /* are the three encoded BCD8 digits followed by a 1-byte length */
1501   /* (significant digits, except that 000 has length 0).  This allows */
1502   /* us to left-align the first declet with non-zero content, then */
1503   /* the remaining ones are full 3-char length.  Fixed-length copies */
1504   /* are used because variable-length memcpy causes a subroutine call */
1505   /* in at least two compilers.  (The copies are length 4 for speed */
1506   /* and are safe because the last item in the array is of length */
1507   /* three and has the length byte following.) */
1508   #define dpd2char(dpdin) u=&DPD2BCD8[((dpdin)&0x3ff)*4];        \
1509          if (c!=cstart) {UINTAT(c)=UINTAT(u)|CHARMASK; c+=3;}    \
1510           else if (*(u+3)) {                                     \
1511            UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; c+=*(u+3);}
1512
1513   #if DECPMAX==7
1514   dpd2char(sourhi>>10);                 /* declet 1 */
1515   dpd2char(sourhi);                     /* declet 2 */
1516
1517   #elif DECPMAX==16
1518   dpd2char(sourhi>>8);                  /* declet 1 */
1519   dpd2char((sourhi<<2) | (sourlo>>30)); /* declet 2 */
1520   dpd2char(sourlo>>20);                 /* declet 3 */
1521   dpd2char(sourlo>>10);                 /* declet 4 */
1522   dpd2char(sourlo);                     /* declet 5 */
1523
1524   #elif DECPMAX==34
1525   dpd2char(sourhi>>4);                  /* declet 1 */
1526   dpd2char((sourhi<<6) | (sourmh>>26)); /* declet 2 */
1527   dpd2char(sourmh>>16);                 /* declet 3 */
1528   dpd2char(sourmh>>6);                  /* declet 4 */
1529   dpd2char((sourmh<<4) | (sourml>>28)); /* declet 5 */
1530   dpd2char(sourml>>18);                 /* declet 6 */
1531   dpd2char(sourml>>8);                  /* declet 7 */
1532   dpd2char((sourml<<2) | (sourlo>>30)); /* declet 8 */
1533   dpd2char(sourlo>>20);                 /* declet 9 */
1534   dpd2char(sourlo>>10);                 /* declet 10 */
1535   dpd2char(sourlo);                     /* declet 11 */
1536   #endif
1537
1538   if (c==cstart) *c++='0';         /* all zeros, empty -- make "0" */
1539
1540   /*[This fast path is valid but adds 3-5 cycles to worst case length] */
1541   /*if (exp==0) {                  // integer or NaN case -- easy */
1542   /*  *c='\0';                     // terminate */
1543   /*  return string; */
1544   /*  } */
1545
1546   e=0;                             /* assume no E */
1547   pre=(Int)(c-cstart)+exp;         /* length+exp  [c->LSD+1] */
1548   /* [here, pre-exp is the digits count (==1 for zero)] */
1549
1550   if (exp>0 || pre<-5) {           /* need exponential form */
1551     e=pre-1;                       /* calculate E value */
1552     pre=1;                         /* assume one digit before '.' */
1553     } /* exponential form */
1554
1555   /* modify the coefficient, adding 0s, '.', and E+nn as needed */
1556   if (pre>0) {                     /* ddd.ddd (plain), perhaps with E */
1557     char *dotat=cstart+pre;
1558     if (dotat<c) {                      /* if embedded dot needed... */
1559       /* move by fours; there must be space for junk at the end */
1560       /* because there is still space for exponent */
1561       s=dotat+ROUNDDOWN4(c-dotat);      /* source */
1562       t=s+1;                            /* target */
1563       /* open the gap */
1564       for (; s>=dotat; s-=4, t-=4) UINTAT(t)=UINTAT(s);
1565       *dotat='.';
1566       c++;                              /* length increased by one */
1567       } /* need dot? */
1568
1569     /* finally add the E-part, if needed; it will never be 0, and has */
1570     /* a maximum length of 3 or 4 digits (asserted above) */
1571     if (e!=0) {
1572       USHORTAT(c)=USHORTAT("E+");       /* starts with E, assume + */
1573       c++;
1574       if (e<0) {
1575         *c='-';                         /* oops, need '-' */
1576         e=-e;                           /* uInt, please */
1577         }
1578       c++;
1579       /* Three-character exponents are easy; 4-character a little trickier */
1580       #if DECEMAXD<=3
1581         u=&BIN2BCD8[e*4];               /* -> 3 digits + length byte */
1582         /* copy fixed 4 characters [is safe], starting at non-zero */
1583         /* and with character mask to convert BCD to char */
1584         UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK;
1585         c+=*(u+3);                      /* bump pointer appropriately */
1586       #elif DECEMAXD==4
1587         if (e<1000) {                   /* 3 (or fewer) digits case */
1588           u=&BIN2BCD8[e*4];             /* -> 3 digits + length byte */
1589           UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; /* [as above] */
1590           c+=*(u+3);                    /* bump pointer appropriately */
1591           }
1592          else {                         /* 4-digits */
1593           Int thou=((e>>3)*1049)>>17;   /* e/1000 */
1594           Int rem=e-(1000*thou);        /* e%1000 */
1595           *c++=(char)('0'+(char)thou);  /* the thousands digit */
1596           u=&BIN2BCD8[rem*4];           /* -> 3 digits + length byte */
1597           UINTAT(c)=UINTAT(u)|CHARMASK; /* copy fixed 3+1 characters [is safe] */
1598           c+=3;                         /* bump pointer, always 3 digits */
1599           }
1600       #endif
1601       }
1602     *c='\0';                            /* add terminator */
1603     /*printf("res %s\n", string); */
1604     return string;
1605     } /* pre>0 */
1606
1607   /* -5<=pre<=0: here for plain 0.ddd or 0.000ddd forms (can never have E) */
1608   /* Surprisingly, this is close to being the worst-case path, so the */
1609   /* shift is done by fours; this is a little tricky because the */
1610   /* rightmost character to be written must not be beyond where the */
1611   /* rightmost terminator could be -- so backoff to not touch */
1612   /* terminator position if need be (this can make exact alignments */
1613   /* for full Doubles, but in some cases needs care not to access too */
1614   /* far to the left) */
1615
1616   pre=-pre+2;                           /* gap width, including "0." */
1617   t=cstart+ROUNDDOWN4(c-cstart)+pre;    /* preferred first target point */
1618   /* backoff if too far to the right */
1619   if (t>string+DECSTRING-5) t=string+DECSTRING-5; /* adjust to fit */
1620   /* now shift the entire coefficient to the right, being careful not */
1621   /* to access to the left of string */
1622   for (s=t-pre; s>=string; s-=4, t-=4) UINTAT(t)=UINTAT(s);
1623   /* for Quads and Singles there may be a character or two left... */
1624   s+=3;                                 /* where next would come from */
1625   for(; s>=cstart; s--, t--) *(t+3)=*(s);
1626   /* now have fill 0. through 0.00000; use overlaps to avoid tests */
1627   if (pre>=4) {
1628     UINTAT(cstart+pre-4)=UINTAT("0000");
1629     UINTAT(cstart)=UINTAT("0.00");
1630     }
1631    else { /* 2 or 3 */
1632     *(cstart+pre-1)='0';
1633     USHORTAT(cstart)=USHORTAT("0.");
1634     }
1635   *(c+pre)='\0';                        /* terminate */
1636   return string;
1637   } /* decFloatToString */
1638
1639 /* ------------------------------------------------------------------ */
1640 /* decFloatToWider -- conversion to next-wider format                 */
1641 /*                                                                    */
1642 /*  source  is the decFloat format number which gets the result of    */
1643 /*          the conversion                                            */
1644 /*  wider   is the decFloatWider format number which will be narrowed */
1645 /*  returns wider                                                     */
1646 /*                                                                    */
1647 /* Widening is always exact; no status is set (sNaNs are copied and   */
1648 /* do not signal).  The result will be canonical if the source is,    */
1649 /* and may or may not be if the source is not.                        */
1650 /* ------------------------------------------------------------------ */
1651 /* widening is not possible for decQuad format numbers; simply omit */
1652 #if !QUAD
1653 decFloatWider * decFloatToWider(const decFloat *source, decFloatWider *wider) {
1654   uInt msd;
1655
1656   /* Construct and copy the sign word */
1657   if (DFISSPECIAL(source)) {
1658     /* copy sign, combination, and first bit of exponent (sNaN selector) */
1659     DFWWORD(wider, 0)=DFWORD(source, 0)&0xfe000000;
1660     msd=0;
1661     }
1662    else { /* is finite number */
1663     uInt exp=GETEXPUN(source)+DECWBIAS; /* get unbiased exponent and rebias */
1664     uInt code=(exp>>DECWECONL)<<29;     /* set two bits of exp [msd=0] */
1665     code|=(exp<<(32-6-DECWECONL)) & 0x03ffffff; /* add exponent continuation */
1666     code|=DFWORD(source, 0)&0x80000000; /* add sign */
1667     DFWWORD(wider, 0)=code;             /* .. and place top word in wider */
1668     msd=GETMSD(source);                 /* get source coefficient MSD [0-9] */
1669     }
1670   /* Copy the coefficient and clear any 'unused' words to left */
1671   #if SINGLE
1672     DFWWORD(wider, 1)=(DFWORD(source, 0)&0x000fffff)|(msd<<20);
1673   #elif DOUBLE
1674     DFWWORD(wider, 2)=(DFWORD(source, 0)&0x0003ffff)|(msd<<18);
1675     DFWWORD(wider, 3)=DFWORD(source, 1);
1676     DFWWORD(wider, 1)=0;
1677   #endif
1678   return wider;
1679   } /* decFloatToWider */
1680 #endif
1681
1682 /* ------------------------------------------------------------------ */
1683 /* decFloatVersion -- return package version string                   */
1684 /*                                                                    */
1685 /*  returns a constant string describing this package                 */
1686 /* ------------------------------------------------------------------ */
1687 const char *decFloatVersion(void) {
1688   return DECVERSION;
1689   } /* decFloatVersion */
1690
1691 /* ------------------------------------------------------------------ */
1692 /* decFloatZero -- set to canonical (integer) zero                    */
1693 /*                                                                    */
1694 /*  df is the decFloat format number to integer +0 (q=0, c=+0)        */
1695 /*  returns df                                                        */
1696 /*                                                                    */
1697 /* No error is possible, and no status can be set.                    */
1698 /* ------------------------------------------------------------------ */
1699 decFloat * decFloatZero(decFloat *df){
1700   DFWORD(df, 0)=ZEROWORD;     /* set appropriate top word */
1701   #if DOUBLE || QUAD
1702     DFWORD(df, 1)=0;
1703     #if QUAD
1704       DFWORD(df, 2)=0;
1705       DFWORD(df, 3)=0;
1706     #endif
1707   #endif
1708   /* decFloatShow(df, "zero"); */
1709   return df;
1710   } /* decFloatZero */
1711
1712 /* ------------------------------------------------------------------ */
1713 /* Private generic function (not format-specific) for development use */
1714 /* ------------------------------------------------------------------ */
1715 /* This is included once only, for all to use */
1716 #if QUAD && (DECCHECK || DECTRACE)
1717   /* ---------------------------------------------------------------- */
1718   /* decShowNum -- display bcd8 number in debug form                  */
1719   /*                                                                  */
1720   /*   num is the bcdnum to display                                   */
1721   /*   tag is a string to label the display                           */
1722   /* ---------------------------------------------------------------- */
1723   void decShowNum(const bcdnum *num, const char *tag) {
1724     const char *csign="+";              /* sign character */
1725     uByte *ub;                          /* work */
1726     if (num->sign==DECFLOAT_Sign) csign="-";
1727
1728     printf(">%s> ", tag);
1729     if (num->exponent==DECFLOAT_Inf) printf("%sInfinity", csign);
1730     else if (num->exponent==DECFLOAT_qNaN) printf("%sqNaN", csign);
1731     else if (num->exponent==DECFLOAT_sNaN) printf("%ssNaN", csign);
1732     else {                              /* finite */
1733      char qbuf[10];                     /* for right-aligned q */
1734      char *c;                           /* work */
1735      const uByte *u;                    /* .. */
1736      Int e=num->exponent;               /* .. exponent */
1737      strcpy(qbuf, "q=");
1738      c=&qbuf[2];                        /* where exponent will go */
1739      /* lay out the exponent */
1740      if (e<0) {
1741        *c++='-';                        /* add '-' */
1742        e=-e;                            /* uInt, please */
1743        }
1744      #if DECEMAXD>4
1745        #error Exponent form is too long for ShowNum to lay out
1746      #endif
1747      if (e==0) *c++='0';                /* 0-length case */
1748       else if (e<1000) {                /* 3 (or fewer) digits case */
1749        u=&BIN2BCD8[e*4];                /* -> 3 digits + length byte */
1750        UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; /* [as above] */
1751        c+=*(u+3);                       /* bump pointer appropriately */
1752        }
1753       else {                            /* 4-digits */
1754        Int thou=((e>>3)*1049)>>17;      /* e/1000 */
1755        Int rem=e-(1000*thou);           /* e%1000 */
1756        *c++=(char)('0'+(char)thou);     /* the thousands digit */
1757        u=&BIN2BCD8[rem*4];              /* -> 3 digits + length byte */
1758        UINTAT(c)=UINTAT(u)|CHARMASK;    /* copy fixed 3+1 characters [is safe] */
1759        c+=3;                            /* bump pointer, always 3 digits */
1760        }
1761      *c='\0';                           /* add terminator */
1762      printf("%7s c=%s", qbuf, csign);
1763      }
1764
1765     if (!EXPISSPECIAL(num->exponent) || num->msd!=num->lsd || *num->lsd!=0) {
1766       for (ub=num->msd; ub<=num->lsd; ub++) { /* coefficient... */
1767         printf("%1x", *ub);
1768         if ((num->lsd-ub)%3==0 && ub!=num->lsd) printf(" "); /* 4-space */
1769         }
1770       }
1771     printf("\n");
1772     } /* decShowNum */
1773 #endif