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