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