decNumber.c, [...]: Upgrade to decNumber 3.61.

	* decNumber.c, decNumber.h, decNumberLocal.h, decDouble.c,
	decDouble.h, decSingle.c, decContext.c, decSingle.h, decPacked.c,
	decCommon.c, decContext.h, decQuad.c, decPacked.h, decQuad.h,
	decDPD.h, decBasic.c: Upgrade to decNumber 3.61.
	* dpd/decimal128.h, dpd/decimal32.c, dpd/decimal32.h,
	dpd/decimal64.c, dpd/decimal128.c, dpd/decimal64.h: Likewise.

From-SVN: r145222

libdecnumber/ChangeLog

+2009-03-29  Ben Elliston  <bje@au.ibm.com>
+
+	* decNumber.c, decNumber.h, decNumberLocal.h, decDouble.c,
+	decDouble.h, decSingle.c, decContext.c, decSingle.h, decPacked.c,
+	decCommon.c, decContext.h, decQuad.c, decPacked.h, decQuad.h,
+	decDPD.h, decBasic.c: Upgrade to decNumber 3.61.
+	* dpd/decimal128.h, dpd/decimal32.c, dpd/decimal32.h,
+	dpd/decimal64.c, dpd/decimal128.c, dpd/decimal64.h: Likewise.
+
 2009-02-10  Joseph Myers  <joseph@codesourcery.com>
 
 	* Makefile.in (clean): Don't remove makedepend$(EXEEXT).

libdecnumber/decBasic.c

@@ -176,7 +176,7 @@ static decFloat * decCanonical(decFloat *result, const decFloat *df) {
 /* op is one of DIVIDE, REMAINDER, DIVIDEINT, or REMNEAR.	      */
 /* ------------------------------------------------------------------ */
 #define DIVCOUNT  0		   /* 1 to instrument subtractions counter */
-#define DIVBASE	  BILLION	   /* the base used for divide */
+#define DIVBASE   ((uInt)BILLION)  /* the base used for divide */
 #define DIVOPLEN  DECPMAX9	   /* operand length ('digits' base 10**9) */
 #define DIVACCLEN (DIVOPLEN*3)	   /* accumulator length (ditto) */
 static decFloat * decDivide(decFloat *result, const decFloat *dfl,
@@ -195,6 +195,7 @@ static decFloat * decDivide(decFloat *result, const decFloat *dfl,
   uInt	 carry, sign;		   /* .. */
   uInt	 *ua, *ud, *uq; 	   /* .. */
   uByte  *ub;			   /* .. */
+  uInt	 uiwork;		   /* for macros */
   uInt	 divtop;		   /* top unit of div adjusted for estimating */
   #if DIVCOUNT
   static uInt maxcount=0;	   /* worst-seen subtractions count */
@@ -421,7 +422,8 @@ static decFloat * decDivide(decFloat *result, const decFloat *dfl,
 	    mul/=*msud;       /* no more to the right */
 	    }
 	   else {
-	    mul=(uLong)(*msua) * (uInt)(DIVBASE<<2) + (*(msua-1)<<2);
+	    mul=(uLong)(*msua) * (uInt)(DIVBASE<<2)
+		+ (*(msua-1)<<2);
 	    mul/=divtop;      /* [divtop already allows for sticky bits] */
 	    }
 	  multiplier=(Int)mul;
@@ -542,8 +544,8 @@ static decFloat * decDivide(decFloat *result, const decFloat *dfl,
   for (uq=msuq, ub=bcdacc+1; uq>=lsuq; uq--, ub+=9) {
     uInt top, mid, rem; 		/* work */
     if (*uq==0) {			/* no split needed */
-      UINTAT(ub)=0;			/* clear 9 BCD8s */
-      UINTAT(ub+4)=0;			/* .. */
+      UBFROMUI(ub, 0);			/* clear 9 BCD8s */
+      UBFROMUI(ub+4, 0);		/* .. */
       *(ub+8)=0;			/* .. */
       continue;
       }
@@ -558,9 +560,9 @@ static decFloat * decDivide(decFloat *result, const decFloat *dfl,
     mid=rem/divsplit6;
     rem=rem%divsplit6;
     /* lay out the nine BCD digits (plus one unwanted byte) */
-    UINTAT(ub)	=UINTAT(&BIN2BCD8[top*4]);
-    UINTAT(ub+3)=UINTAT(&BIN2BCD8[mid*4]);
-    UINTAT(ub+6)=UINTAT(&BIN2BCD8[rem*4]);
+    UBFROMUI(ub,   UBTOUI(&BIN2BCD8[top*4]));
+    UBFROMUI(ub+3, UBTOUI(&BIN2BCD8[mid*4]));
+    UBFROMUI(ub+6, UBTOUI(&BIN2BCD8[rem*4]));
     } /* BCD conversion loop */
   ub--; 				/* -> lsu */
 
@@ -651,7 +653,7 @@ static decFloat * decDivide(decFloat *result, const decFloat *dfl,
 	if (bump!=0) {			     /* need increment */
 	  /* increment the coefficient; this might end up with 1000... */
 	  ub=num.lsd;
-	  for (; UINTAT(ub-3)==0x09090909; ub-=4) UINTAT(ub-3)=0;
+	  for (; UBTOUI(ub-3)==0x09090909; ub-=4) UBFROMUI(ub-3, 0);
 	  for (; *ub==9; ub--) *ub=0;	     /* at most 3 more */
 	  *ub+=1;
 	  if (ub<num.msd) num.msd--;	     /* carried */
@@ -704,7 +706,7 @@ static decFloat * decDivide(decFloat *result, const decFloat *dfl,
 /* during lazy carry splitting because the initial quotient estimate */
 /* (est) can exceed 32 bits. */
 
-#define MULTBASE  BILLION	   /* the base used for multiply */
+#define MULTBASE  ((uInt)BILLION)  /* the base used for multiply */
 #define MULOPLEN  DECPMAX9	   /* operand length ('digits' base 10**9) */
 #define MULACCLEN (MULOPLEN*2)		    /* accumulator length (ditto) */
 #define LEADZEROS (MULACCLEN*9 - DECPMAX*2) /* leading zeros always */
@@ -724,6 +726,7 @@ static void decFiniteMultiply(bcdnum *num, uByte *bcdacc,
   uInt	 bufr[MULOPLEN];	   /* right coefficient (base-billion) */
   uInt	 *ui, *uj;		   /* work */
   uByte  *ub;			   /* .. */
+  uInt	 uiwork;		   /* for macros */
 
   #if DECUSE64
   uLong  accl[MULACCLEN];	   /* lazy accumulator (base-billion+) */
@@ -840,7 +843,7 @@ static void decFiniteMultiply(bcdnum *num, uByte *bcdacc,
   for (pl=accl, pa=acc; pl<accl+MULACCLEN; pl++, pa++) { /* each column position */
     uInt lo, hop;			/* work */
     uInt est;				/* cannot exceed 4E+9 */
-    if (*pl>MULTBASE) {
+    if (*pl>=MULTBASE) {
       /* *pl holds a binary number which needs to be split */
       hop=(uInt)(*pl>>MULSHIFTA);
       est=(uInt)(((uLong)hop*MULMAGIC)>>MULSHIFTB);
@@ -1042,13 +1045,13 @@ static void decFiniteMultiply(bcdnum *num, uByte *bcdacc,
       mid=rem/mulsplit6;
       rem=rem%mulsplit6;
       /* lay out the nine BCD digits (plus one unwanted byte) */
-      UINTAT(ub)  =UINTAT(&BIN2BCD8[top*4]);
-      UINTAT(ub+3)=UINTAT(&BIN2BCD8[mid*4]);
-      UINTAT(ub+6)=UINTAT(&BIN2BCD8[rem*4]);
+      UBFROMUI(ub,   UBTOUI(&BIN2BCD8[top*4]));
+      UBFROMUI(ub+3, UBTOUI(&BIN2BCD8[mid*4]));
+      UBFROMUI(ub+6, UBTOUI(&BIN2BCD8[rem*4]));
       }
      else {				/* *pa==0 */
-      UINTAT(ub)=0;			/* clear 9 BCD8s */
-      UINTAT(ub+4)=0;			/* .. */
+      UBFROMUI(ub, 0);			/* clear 9 BCD8s */
+      UBFROMUI(ub+4, 0);		/* .. */
       *(ub+8)=0;			/* .. */
       }
     if (pa==acc) break;
@@ -1096,20 +1099,39 @@ decFloat * decFloatAbs(decFloat *result, const decFloat *df,
 /*   returns result						      */
 /*								      */
 /* ------------------------------------------------------------------ */
+#if QUAD
+/* Table for testing MSDs for fastpath elimination; returns the MSD of */
+/* a decDouble or decQuad (top 6 bits tested) ignoring the sign. */
+/* Infinities return -32 and NaNs return -128 so that summing the two */
+/* MSDs also allows rapid tests for the Specials (see code below). */
+const Int DECTESTMSD[64]={
+  0, 1, 2, 3, 4, 5, 6, 7, 0, 1, 2, 3, 4, 5,   6,    7,
+  0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 9, 8, 9, -32, -128,
+  0, 1, 2, 3, 4, 5, 6, 7, 0, 1, 2, 3, 4, 5,   6,    7,
+  0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 9, 8, 9, -32, -128};
+#else
+/* The table for testing MSDs is shared between the modules */
+extern const Int DECTESTMSD[64];
+#endif
+
 decFloat * decFloatAdd(decFloat *result,
 		       const decFloat *dfl, const decFloat *dfr,
 		       decContext *set) {
   bcdnum num;			   /* for final conversion */
-  Int	 expl, expr;		   /* left and right exponents */
-  uInt	 *ui, *uj;		   /* work */
-  uByte	 *ub;			   /* .. */
+  Int	 bexpl, bexpr;		   /* left and right biased exponents */
+  uByte  *ub, *us, *ut; 	   /* work */
+  uInt	 uiwork;		   /* for macros */
+  #if QUAD
+  uShort uswork;		   /* .. */
+  #endif
 
   uInt sourhil, sourhir;	   /* top words from source decFloats */
-				   /* [valid only until specials */
-				   /* handled or exponents decoded] */
+				   /* [valid only through end of */
+				   /* fastpath code -- before swap] */
   uInt diffsign;		   /* non-zero if signs differ */
   uInt carry;			   /* carry: 0 or 1 before add loop */
   Int  overlap; 		   /* coefficient overlap (if full) */
+  Int  summ;			   /* sum of the MSDs */
   /* the following buffers hold coefficients with various alignments */
   /* (see commentary and diagrams below) */
   uByte acc[4+2+DECPMAX*3+8];
@@ -1123,42 +1145,110 @@ decFloat * decFloatAdd(decFloat *result,
   #endif
 
   /* Start decoding the arguments */
-  /* the initial exponents are placed into the opposite Ints to */
+  /* The initial exponents are placed into the opposite Ints to */
   /* that which might be expected; there are two sets of data to */
   /* keep track of (each decFloat and the corresponding exponent), */
   /* and this scheme means that at the swap point (after comparing */
   /* exponents) only one pair of words needs to be swapped */
-  /* whichever path is taken (thereby minimising worst-case path) */
+  /* whichever path is taken (thereby minimising worst-case path). */
+  /* The calculated exponents will be nonsense when the arguments are */
+  /* Special, but are not used in that path */
   sourhil=DFWORD(dfl, 0);	   /* LHS top word */
-  expr=DECCOMBEXP[sourhil>>26];	   /* get exponent high bits (in place) */
+  summ=DECTESTMSD[sourhil>>26];    /* get first MSD for testing */
+  bexpr=DECCOMBEXP[sourhil>>26];   /* get exponent high bits (in place) */
+  bexpr+=GETECON(dfl);		   /* .. + continuation */
+
   sourhir=DFWORD(dfr, 0);	   /* RHS top word */
-  expl=DECCOMBEXP[sourhir>>26];
+  summ+=DECTESTMSD[sourhir>>26];   /* sum MSDs for testing */
+  bexpl=DECCOMBEXP[sourhir>>26];
+  bexpl+=GETECON(dfr);
+
+  /* here bexpr has biased exponent from lhs, and vice versa */
 
   diffsign=(sourhil^sourhir)&DECFLOAT_Sign;
 
-  if (EXPISSPECIAL(expl | expr)) { /* either is special? */
-    if (DFISNAN(dfl) || DFISNAN(dfr)) return decNaNs(result, dfl, dfr, set);
-    /* one or two infinities */
+  /* now determine whether to take a fast path or the full-function */
+  /* slow path.  The slow path must be taken when: */
+  /*   -- both numbers are finite, and: */
+  /*	     the exponents are different, or */
+  /*	     the signs are different, or */
+  /*	     the sum of the MSDs is >8 (hence might overflow) */
+  /* specialness and the sum of the MSDs can be tested at once using */
+  /* the summ value just calculated, so the test for specials is no */
+  /* longer on the worst-case path (as of 3.60) */
+
+  if (summ<=8) {		   /* MSD+MSD is good, or there is a special */
+    if (summ<0) {		   /* there is a special */
+      /* Inf+Inf would give -64; Inf+finite is -32 or higher */
+      if (summ<-64) return decNaNs(result, dfl, dfr, set);  /* one or two NaNs */
       /* two infinities with different signs is invalid */
-    if (diffsign && DFISINF(dfl) && DFISINF(dfr))
-      return decInvalid(result, set);
+      if (summ==-64 && diffsign) return decInvalid(result, set);
       if (DFISINF(dfl)) return decInfinity(result, dfl);    /* LHS is infinite */
-    return decInfinity(result, dfr);		       /* RHS must be Infinite */
+      return decInfinity(result, dfr);			    /* RHS must be Inf */
       }
+    /* Here when both arguments are finite; fast path is possible */
+    /* (currently only for aligned and same-sign) */
+    if (bexpr==bexpl && !diffsign) {
+      uInt tac[DECLETS+1];		/* base-1000 coefficient */
+      uInt encode;			/* work */
 
-  /* Here when both arguments are finite */
+      /* Get one coefficient as base-1000 and add the other */
+      GETCOEFFTHOU(dfl, tac);		/* least-significant goes to [0] */
+      ADDCOEFFTHOU(dfr, tac);
+      /* here the sum of the MSDs (plus any carry) will be <10 due to */
+      /* the fastpath test earlier */
+
+      /* construct the result; low word is the same for both formats */
+      encode =BIN2DPD[tac[0]];
+      encode|=BIN2DPD[tac[1]]<<10;
+      encode|=BIN2DPD[tac[2]]<<20;
+      encode|=BIN2DPD[tac[3]]<<30;
+      DFWORD(result, (DECBYTES/4)-1)=encode;
+
+      /* collect next two declets (all that remains, for Double) */
+      encode =BIN2DPD[tac[3]]>>2;
+      encode|=BIN2DPD[tac[4]]<<8;
+
+      #if QUAD
+      /* complete and lay out middling words */
+      encode|=BIN2DPD[tac[5]]<<18;
+      encode|=BIN2DPD[tac[6]]<<28;
+      DFWORD(result, 2)=encode;
 
-  /* complete exponent gathering (keeping swapped) */
-  expr+=GETECON(dfl)-DECBIAS;	   /* .. + continuation and unbias */
-  expl+=GETECON(dfr)-DECBIAS;
-  /* here expr has exponent from lhs, and vice versa */
+      encode =BIN2DPD[tac[6]]>>4;
+      encode|=BIN2DPD[tac[7]]<<6;
+      encode|=BIN2DPD[tac[8]]<<16;
+      encode|=BIN2DPD[tac[9]]<<26;
+      DFWORD(result, 1)=encode;
+
+      /* and final two declets */
+      encode =BIN2DPD[tac[9]]>>6;
+      encode|=BIN2DPD[tac[10]]<<4;
+      #endif
+
+      /* add exponent continuation and sign (from either argument) */
+      encode|=sourhil & (ECONMASK | DECFLOAT_Sign);
+
+      /* create lookup index = MSD + top two bits of biased exponent <<4 */
+      tac[DECLETS]|=(bexpl>>DECECONL)<<4;
+      encode|=DECCOMBFROM[tac[DECLETS]]; /* add constructed combination field */
+      DFWORD(result, 0)=encode; 	 /* complete */
+
+      /* decFloatShow(result, ">"); */
+      return result;
+      } /* fast path OK */
+    /* drop through to slow path */
+    } /* low sum or Special(s) */
+
+  /* Slow path required -- arguments are finite and might overflow,   */
+  /* or require alignment, or might have different signs	      */
 
   /* now swap either exponents or argument pointers */
-  if (expl<=expr) {
+  if (bexpl<=bexpr) {
     /* original left is bigger */
-    Int expswap=expl;
-    expl=expr;
-    expr=expswap;
+    Int bexpswap=bexpl;
+    bexpl=bexpr;
+    bexpr=bexpswap;
     /* printf("left bigger\n"); */
     }
    else {
@@ -1167,7 +1257,7 @@ decFloat * decFloatAdd(decFloat *result,
     dfr=dfswap;
     /* printf("right bigger\n"); */
     }
-  /* [here dfl and expl refer to the datum with the larger exponent, */
+  /* [here dfl and bexpl refer to the datum with the larger exponent, */
   /* of if the exponents are equal then the original LHS argument] */
 
   /* if lhs is zero then result will be the rhs (now known to have */
@@ -1209,19 +1299,19 @@ decFloat * decFloatAdd(decFloat *result,
   #if DOUBLE
     #define COFF 4			/* offset into acc */
   #elif QUAD
-    USHORTAT(acc+4)=0;			/* prefix 00 */
+    UBFROMUS(acc+4, 0); 		/* prefix 00 */
     #define COFF 6			/* offset into acc */
   #endif
 
   GETCOEFF(dfl, acc+COFF);		/* decode from decFloat */
   ulsd=acc+COFF+DECPMAX-1;
   umsd=acc+4;				/* [having this here avoids */
-					/* weird GCC optimizer failure] */
+
   #if DECTRACE
   {bcdnum tum;
   tum.msd=umsd;
   tum.lsd=ulsd;
-  tum.exponent=expl;
+  tum.exponent=bexpl-DECBIAS;
   tum.sign=DFWORD(dfl, 0) & DECFLOAT_Sign;
   decShowNum(&tum, "dflx");}
   #endif
@@ -1235,16 +1325,16 @@ decFloat * decFloatAdd(decFloat *result,
   carry=0;				/* assume no carry */
   if (diffsign) {
     carry=CARRYPAT;			/* for +1 during add */
-    UINTAT(acc+ 4)=0x09090909-UINTAT(acc+ 4);
-    UINTAT(acc+ 8)=0x09090909-UINTAT(acc+ 8);
-    UINTAT(acc+12)=0x09090909-UINTAT(acc+12);
-    UINTAT(acc+16)=0x09090909-UINTAT(acc+16);
+    UBFROMUI(acc+ 4, 0x09090909-UBTOUI(acc+ 4));
+    UBFROMUI(acc+ 8, 0x09090909-UBTOUI(acc+ 8));
+    UBFROMUI(acc+12, 0x09090909-UBTOUI(acc+12));
+    UBFROMUI(acc+16, 0x09090909-UBTOUI(acc+16));
     #if QUAD
-    UINTAT(acc+20)=0x09090909-UINTAT(acc+20);
-    UINTAT(acc+24)=0x09090909-UINTAT(acc+24);
-    UINTAT(acc+28)=0x09090909-UINTAT(acc+28);
-    UINTAT(acc+32)=0x09090909-UINTAT(acc+32);
-    UINTAT(acc+36)=0x09090909-UINTAT(acc+36);
+    UBFROMUI(acc+20, 0x09090909-UBTOUI(acc+20));
+    UBFROMUI(acc+24, 0x09090909-UBTOUI(acc+24));
+    UBFROMUI(acc+28, 0x09090909-UBTOUI(acc+28));
+    UBFROMUI(acc+32, 0x09090909-UBTOUI(acc+32));
+    UBFROMUI(acc+36, 0x09090909-UBTOUI(acc+36));
     #endif
     } /* diffsign */
 
@@ -1252,9 +1342,9 @@ decFloat * decFloatAdd(decFloat *result,
   /* it can be put straight into acc (with an appropriate gap, if */
   /* needed) because no actual addition will be needed (except */
   /* possibly to complete ten's complement) */
-  overlap=DECPMAX-(expl-expr);
+  overlap=DECPMAX-(bexpl-bexpr);
   #if DECTRACE
-  printf("exps: %ld %ld\n", (LI)expl, (LI)expr);
+  printf("exps: %ld %ld\n", (LI)(bexpl-DECBIAS), (LI)(bexpr-DECBIAS));
   printf("Overlap=%ld carry=%08lx\n", (LI)overlap, (LI)carry);
   #endif
 
@@ -1274,13 +1364,13 @@ decFloat * decFloatAdd(decFloat *result,
     /* safe because the lhs is non-zero]. */
     gap=-overlap;
     if (gap>DECPMAX) {
-      expr+=gap-1;
+      bexpr+=gap-1;
       gap=DECPMAX;
       }
     ub=ulsd+gap+1;			/* where MSD will go */
     /* Fill the gap with 0s; note that there is no addition to do */
-    ui=&UINTAT(acc+COFF+DECPMAX);	/* start of gap */
-    for (; ui<&UINTAT(ub); ui++) *ui=0; /* mind the gap */
+    ut=acc+COFF+DECPMAX;		/* start of gap */
+    for (; ut<ub; ut+=4) UBFROMUI(ut, 0); /* mind the gap */
     if (overlap<-DECPMAX) {		/* gap was > DECPMAX */
       *ub=(uByte)(!DFISZERO(dfr));	/* make sticky digit */
       }
@@ -1294,63 +1384,74 @@ decFloat * decFloatAdd(decFloat *result,
    else {				/* overlap>0 */
     /* coefficients overlap (perhaps completely, although also */
     /* perhaps only where zeros) */
+    if (overlap==DECPMAX) {		/* aligned */
+      ub=buf+COFF;			/* where msd will go */
+      #if QUAD
+      UBFROMUS(buf+4, 0);		/* clear quad's 00 */
+      #endif
+      GETCOEFF(dfr, ub);		/* decode from decFloat */
+      }
+     else {				/* unaligned */
       ub=buf+COFF+DECPMAX-overlap;	/* where MSD will go */
       /* Fill the prefix gap with 0s; 8 will cover most common */
       /* unalignments, so start with direct assignments (a loop is */
       /* then used for any remaining -- the loop (and the one in a */
       /* moment) is not then on the critical path because the number */
       /* of additions is reduced by (at least) two in this case) */
-    UINTAT(buf+4)=0;			/* [clears decQuad 00 too] */
-    UINTAT(buf+8)=0;
+      UBFROMUI(buf+4, 0);		/* [clears decQuad 00 too] */
+      UBFROMUI(buf+8, 0);
       if (ub>buf+12) {
-      ui=&UINTAT(buf+12);		/* start of any remaining */
-      for (; ui<&UINTAT(ub); ui++) *ui=0; /* fill them */
+	ut=buf+12;			/* start any remaining */
+	for (; ut<ub; ut+=4) UBFROMUI(ut, 0); /* fill them */
 	}
       GETCOEFF(dfr, ub);		/* decode from decFloat */
 
       /* now move tail of rhs across to main acc; again use direct */
-    /* assignment for 8 digits-worth */
-    UINTAT(acc+COFF+DECPMAX)=UINTAT(buf+COFF+DECPMAX);
-    UINTAT(acc+COFF+DECPMAX+4)=UINTAT(buf+COFF+DECPMAX+4);
+      /* copies for 8 digits-worth */
+      UBFROMUI(acc+COFF+DECPMAX,   UBTOUI(buf+COFF+DECPMAX));
+      UBFROMUI(acc+COFF+DECPMAX+4, UBTOUI(buf+COFF+DECPMAX+4));
       if (buf+COFF+DECPMAX+8<ub+DECPMAX) {
-      uj=&UINTAT(buf+COFF+DECPMAX+8);	/* source */
-      ui=&UINTAT(acc+COFF+DECPMAX+8);	/* target */
-      for (; uj<&UINTAT(ub+DECPMAX); ui++, uj++) *ui=*uj;
+	us=buf+COFF+DECPMAX+8;		/* source */
+	ut=acc+COFF+DECPMAX+8;		/* target */
+	for (; us<ub+DECPMAX; us+=4, ut+=4) UBFROMUI(ut, UBTOUI(us));
 	}
+      } /* unaligned */
 
     ulsd=acc+(ub-buf+DECPMAX-1);	/* update LSD pointer */
 
-    /* now do the add of the non-tail; this is all nicely aligned, */
+    /* Now do the add of the non-tail; this is all nicely aligned, */
     /* and is over a multiple of four digits (because for Quad two */
-    /* two 0 digits were added on the left); words in both acc and */
+    /* zero digits were added on the left); words in both acc and */
     /* buf (buf especially) will often be zero */
-    /* [byte-by-byte add, here, is about 15% slower than the by-fours] */
+    /* [byte-by-byte add, here, is about 15% slower total effect than */
+    /* the by-fours] */
 
     /* Now effect the add; this is harder on a little-endian */
     /* machine as the inter-digit carry cannot use the usual BCD */
     /* addition trick because the bytes are loaded in the wrong order */
     /* [this loop could be unrolled, but probably scarcely worth it] */
 
-    ui=&UINTAT(acc+COFF+DECPMAX-4);	/* target LSW (acc) */
-    uj=&UINTAT(buf+COFF+DECPMAX-4);	/* source LSW (buf, to add to acc) */
+    ut=acc+COFF+DECPMAX-4;		/* target LSW (acc) */
+    us=buf+COFF+DECPMAX-4;		/* source LSW (buf, to add to acc) */
 
     #if !DECLITEND
-    for (; ui>=&UINTAT(acc+4); ui--, uj--) {
+    for (; ut>=acc+4; ut-=4, us-=4) {	/* big-endian add loop */
       /* bcd8 add */
-      carry+=*uj;			/* rhs + carry */
+      carry+=UBTOUI(us);		/* rhs + carry */
       if (carry==0) continue;		/* no-op */
-      carry+=*ui;			/* lhs */
+      carry+=UBTOUI(ut);		/* lhs */
       /* Big-endian BCD adjust (uses internal carry) */
       carry+=0x76f6f6f6;		/* note top nibble not all bits */
-      *ui=(carry & 0x0f0f0f0f) - ((carry & 0x60606060)>>4); /* BCD adjust */
+      /* apply BCD adjust and save */
+      UBFROMUI(ut, (carry & 0x0f0f0f0f) - ((carry & 0x60606060)>>4));
       carry>>=31;			/* true carry was at far left */
       } /* add loop */
     #else
-    for (; ui>=&UINTAT(acc+4); ui--, uj--) {
+    for (; ut>=acc+4; ut-=4, us-=4) {	/* little-endian add loop */
       /* bcd8 add */
-      carry+=*uj;			/* rhs + carry */
+      carry+=UBTOUI(us);		/* rhs + carry */
       if (carry==0) continue;		/* no-op [common if unaligned] */
-      carry+=*ui;			/* lhs */
+      carry+=UBTOUI(ut);		/* lhs */
       /* Little-endian BCD adjust; inter-digit carry must be manual */
       /* because the lsb from the array will be in the most-significant */
       /* byte of carry */
@@ -1359,12 +1460,13 @@ decFloat * decFloatAdd(decFloat *result,
       carry+=(carry & 0x00800000)>>15;
       carry+=(carry & 0x00008000)>>15;
       carry-=(carry & 0x60606060)>>4;	/* BCD adjust back */
-      *ui=carry & 0x0f0f0f0f;		/* clear debris and save */
+      UBFROMUI(ut, carry & 0x0f0f0f0f); /* clear debris and save */
       /* here, final carry-out bit is at 0x00000080; move it ready */
       /* for next word-add (i.e., to 0x01000000) */
       carry=(carry & 0x00000080)<<17;
       } /* add loop */
     #endif
+
     #if DECTRACE
     {bcdnum tum;
     printf("Add done, carry=%08lx, diffsign=%ld\n", (LI)carry, (LI)diffsign);
@@ -1392,36 +1494,36 @@ decFloat * decFloatAdd(decFloat *result,
       *(ulsd+1)=0;
       #endif
       /* there are always at least four coefficient words */
-      UINTAT(umsd)   =0x09090909-UINTAT(umsd);
-      UINTAT(umsd+4) =0x09090909-UINTAT(umsd+4);
-      UINTAT(umsd+8) =0x09090909-UINTAT(umsd+8);
-      UINTAT(umsd+12)=0x09090909-UINTAT(umsd+12);
+      UBFROMUI(umsd,	0x09090909-UBTOUI(umsd));
+      UBFROMUI(umsd+4,	0x09090909-UBTOUI(umsd+4));
+      UBFROMUI(umsd+8,	0x09090909-UBTOUI(umsd+8));
+      UBFROMUI(umsd+12, 0x09090909-UBTOUI(umsd+12));
       #if DOUBLE
 	#define BNEXT 16
       #elif QUAD
-	UINTAT(umsd+16)=0x09090909-UINTAT(umsd+16);
-	UINTAT(umsd+20)=0x09090909-UINTAT(umsd+20);
-	UINTAT(umsd+24)=0x09090909-UINTAT(umsd+24);
-	UINTAT(umsd+28)=0x09090909-UINTAT(umsd+28);
-	UINTAT(umsd+32)=0x09090909-UINTAT(umsd+32);
+	UBFROMUI(umsd+16, 0x09090909-UBTOUI(umsd+16));
+	UBFROMUI(umsd+20, 0x09090909-UBTOUI(umsd+20));
+	UBFROMUI(umsd+24, 0x09090909-UBTOUI(umsd+24));
+	UBFROMUI(umsd+28, 0x09090909-UBTOUI(umsd+28));
+	UBFROMUI(umsd+32, 0x09090909-UBTOUI(umsd+32));
 	#define BNEXT 36
       #endif
       if (ulsd>=umsd+BNEXT) {		/* unaligned */
 	/* eight will handle most unaligments for Double; 16 for Quad */
-	UINTAT(umsd+BNEXT)=0x09090909-UINTAT(umsd+BNEXT);
-	UINTAT(umsd+BNEXT+4)=0x09090909-UINTAT(umsd+BNEXT+4);
+	UBFROMUI(umsd+BNEXT,   0x09090909-UBTOUI(umsd+BNEXT));
+	UBFROMUI(umsd+BNEXT+4, 0x09090909-UBTOUI(umsd+BNEXT+4));
 	#if DOUBLE
 	#define BNEXTY (BNEXT+8)
 	#elif QUAD
-	UINTAT(umsd+BNEXT+8)=0x09090909-UINTAT(umsd+BNEXT+8);
-	UINTAT(umsd+BNEXT+12)=0x09090909-UINTAT(umsd+BNEXT+12);
+	UBFROMUI(umsd+BNEXT+8,	0x09090909-UBTOUI(umsd+BNEXT+8));
+	UBFROMUI(umsd+BNEXT+12, 0x09090909-UBTOUI(umsd+BNEXT+12));
 	#define BNEXTY (BNEXT+16)
 	#endif
 	if (ulsd>=umsd+BNEXTY) {	/* very unaligned */
-	  ui=&UINTAT(umsd+BNEXTY);	/* -> continue */
-	  for (;;ui++) {
-	    *ui=0x09090909-*ui;		/* invert four digits */
-	    if (ui>=&UINTAT(ulsd-3)) break; /* all done */
+	  ut=umsd+BNEXTY;		/* -> continue */
+	  for (;;ut+=4) {
+	    UBFROMUI(ut, 0x09090909-UBTOUI(ut)); /* invert four digits */
+	    if (ut>=ulsd-3) break;	/* all done */
 	    }
 	  }
 	}
@@ -1446,7 +1548,7 @@ decFloat * decFloatAdd(decFloat *result,
 	umsd=acc+COFF+DECPMAX-1;   /* so far, so zero */
 	if (ulsd>umsd) {	   /* more to check */
 	  umsd++;		   /* to align after checked area */
-	  for (; UINTAT(umsd)==0 && umsd+3<ulsd;) umsd+=4;
+	  for (; UBTOUI(umsd)==0 && umsd+3<ulsd;) umsd+=4;
 	  for (; *umsd==0 && umsd<ulsd;) umsd++;
 	  }
 	if (*umsd==0) { 	   /* must be true zero (and diffsign) */
@@ -1470,7 +1572,7 @@ decFloat * decFloatAdd(decFloat *result,
 
   num.msd=umsd; 		   /* set MSD .. */
   num.lsd=ulsd; 		   /* .. and LSD */
-  num.exponent=expr;		   /* set exponent to smaller */
+  num.exponent=bexpr-DECBIAS;	   /* set exponent to smaller, unbiassed */
 
   #if DECTRACE
   decFloatShow(dfl, "dfl");
@@ -1640,6 +1742,10 @@ decFloat * decFloatCompareSignal(decFloat *result,
 decFloat * decFloatCompareTotal(decFloat *result,
 				const decFloat *dfl, const decFloat *dfr) {
   Int  comp;				     /* work */
+  uInt uiwork;				     /* for macros */
+  #if QUAD
+  uShort uswork;			     /* .. */
+  #endif
   if (DFISNAN(dfl) || DFISNAN(dfr)) {
     Int nanl, nanr;			     /* work */
     /* morph NaNs to +/- 1 or 2, leave numbers as 0 */
@@ -1660,17 +1766,16 @@ decFloat * decFloatCompareTotal(decFloat *result,
       /* decode the coefficients */
       /* (shift both right two if Quad to make a multiple of four) */
       #if QUAD
-	ub = bufl;                           /* avoid type-pun violation */
-	USHORTAT(ub)=0;
-	uc = bufr;                           /* avoid type-pun violation */
-	USHORTAT(uc)=0;
+	UBFROMUS(bufl, 0);
+	UBFROMUS(bufr, 0);
       #endif
       GETCOEFF(dfl, bufl+QUAD*2);	     /* decode from decFloat */
       GETCOEFF(dfr, bufr+QUAD*2);	     /* .. */
       /* all multiples of four, here */
       comp=0;				     /* assume equal */
       for (ub=bufl, uc=bufr; ub<bufl+DECPMAX+QUAD*2; ub+=4, uc+=4) {
-	if (UINTAT(ub)==UINTAT(uc)) continue; /* so far so same */
+	uInt ui=UBTOUI(ub);
+	if (ui==UBTOUI(uc)) continue; /* so far so same */
 	/* about to find a winner; go by bytes in case little-endian */
 	for (;; ub++, uc++) {
 	  if (*ub==*uc) continue;
@@ -1906,21 +2011,23 @@ decFloat * decFloatDivideInteger(decFloat *result,
 decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
 		       const decFloat *dfr, const decFloat *dff,
 		       decContext *set) {
+
   /* The accumulator has the bytes needed for FiniteMultiply, plus */
   /* one byte to the left in case of carry, plus DECPMAX+2 to the */
   /* right for the final addition (up to full fhs + round & sticky) */
-  #define FMALEN (1+ (DECPMAX9*18) +DECPMAX+2)
+  #define FMALEN (ROUNDUP4(1+ (DECPMAX9*18+1) +DECPMAX+2))
   uByte  acc[FMALEN];		   /* for multiplied coefficient in BCD */
 				   /* .. and for final result */
   bcdnum mul;			   /* for multiplication result */
   bcdnum fin;			   /* for final operand, expanded */
-  uByte	 coe[DECPMAX];		   /* dff coefficient in BCD */
+  uByte  coe[ROUNDUP4(DECPMAX)];   /* dff coefficient in BCD */
   bcdnum *hi, *lo;		   /* bcdnum with higher/lower exponent */
   uInt	 diffsign;		   /* non-zero if signs differ */
   uInt	 hipad; 		   /* pad digit for hi if needed */
   Int	 padding;		   /* excess exponent */
   uInt	 carry; 		   /* +1 for ten's complement and during add */
   uByte  *ub, *uh, *ul; 	   /* work */
+  uInt	 uiwork;		   /* for macros */
 
   /* handle all the special values [any special operand leads to a */
   /* special result] */
@@ -1971,8 +2078,8 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
   GETCOEFF(dff, coe);			/* extract the coefficient */
 
   /* now set hi and lo so that hi points to whichever of mul and fin */
-  /* has the higher exponent and lo point to the other [don't care if */
-  /* the same] */
+  /* has the higher exponent and lo points to the other [don't care, */
+  /* if the same].  One coefficient will be in acc, the other in coe. */
   if (mul.exponent>=fin.exponent) {
     hi=&mul;
     lo=&fin;
@@ -1983,22 +2090,23 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
     }
 
   /* remove leading zeros on both operands; this will save time later */
-  /* and make testing for zero trivial */
-  for (; UINTAT(hi->msd)==0 && hi->msd+3<hi->lsd;) hi->msd+=4;
+  /* and make testing for zero trivial (tests are safe because acc */
+  /* and coe are rounded up to uInts) */
+  for (; UBTOUI(hi->msd)==0 && hi->msd+3<hi->lsd;) hi->msd+=4;
   for (; *hi->msd==0 && hi->msd<hi->lsd;) hi->msd++;
-  for (; UINTAT(lo->msd)==0 && lo->msd+3<lo->lsd;) lo->msd+=4;
+  for (; UBTOUI(lo->msd)==0 && lo->msd+3<lo->lsd;) lo->msd+=4;
   for (; *lo->msd==0 && lo->msd<lo->lsd;) lo->msd++;
 
   /* if hi is zero then result will be lo (which has the smaller */
   /* exponent), which also may need to be tested for zero for the */
   /* weird IEEE 754 sign rules */
-  if (*hi->msd==0 && hi->msd==hi->lsd) {     /* hi is zero */
+  if (*hi->msd==0) {			     /* hi is zero */
     /* "When the sum of two operands with opposite signs is */
     /* exactly zero, the sign of that sum shall be '+' in all */
     /* rounding modes except round toward -Infinity, in which */
     /* mode that sign shall be '-'." */
     if (diffsign) {
-      if (*lo->msd==0 && lo->msd==lo->lsd) { /* lo is zero */
+      if (*lo->msd==0) {		     /* lo is zero */
 	lo->sign=0;
 	if (set->round==DEC_ROUND_FLOOR) lo->sign=DECFLOAT_Sign;
 	} /* diffsign && lo=0 */
@@ -2006,10 +2114,11 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
     return decFinalize(result, lo, set);     /* may need clamping */
     } /* numfl is zero */
   /* [here, both are minimal length and hi is non-zero] */
+  /* (if lo is zero then padding with zeros may be needed, below) */
 
   /* if signs differ, take the ten's complement of hi (zeros to the */
-  /* right do not matter because the complement of zero is zero); */
-  /* the +1 is done later, as part of the addition, inserted at the */
+  /* right do not matter because the complement of zero is zero); the */
+  /* +1 is done later, as part of the addition, inserted at the */
   /* correct digit */
   hipad=0;
   carry=0;
@@ -2017,7 +2126,7 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
     hipad=9;
     carry=1;
     /* exactly the correct number of digits must be inverted */
-    for (uh=hi->msd; uh<hi->lsd-3; uh+=4) UINTAT(uh)=0x09090909-UINTAT(uh);
+    for (uh=hi->msd; uh<hi->lsd-3; uh+=4) UBFROMUI(uh, 0x09090909-UBTOUI(uh));
     for (; uh<=hi->lsd; uh++) *uh=(uByte)(0x09-*uh);
     }
 
@@ -2032,7 +2141,8 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
   /* printf("FMA pad %ld\n", (LI)padding); */
 
   /* the result of the addition will be built into the accumulator, */
-  /* starting from the far right; this could be either hi or lo */
+  /* starting from the far right; this could be either hi or lo, and */
+  /* will be aligned */
   ub=acc+FMALEN-1;		   /* where lsd of result will go */
   ul=lo->lsd;			   /* lsd of rhs */
 
@@ -2042,45 +2152,43 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
     /* digit at the right place, as it stays clear of hi digits */
     /* [it must be DECPMAX+2 because during a subtraction the msd */
     /* could become 0 after a borrow from 1.000 to 0.9999...] */
-    Int hilen=(Int)(hi->lsd-hi->msd+1); /* lengths */
-    Int lolen=(Int)(lo->lsd-lo->msd+1); /* .. */
+
+    Int hilen=(Int)(hi->lsd-hi->msd+1); /* length of hi */
+    Int lolen=(Int)(lo->lsd-lo->msd+1); /* and of lo */
+
+    if (hilen+padding-lolen > DECPMAX+2) {   /* can reduce lo to single */
+      /* make sure it is virtually at least DECPMAX from hi->msd, at */
+      /* least to right of hi->lsd (in case of destructive subtract), */
+      /* and separated by at least two digits from either of those */
+      /* (the tricky DOUBLE case is when hi is a 1 that will become a */
+      /* 0.9999... by subtraction: */
+      /*   hi:	 1				     E+16 */
+      /*   lo:	  .................1000000000000000  E-16 */
+      /* which for the addition pads to: */
+      /*   hi:	 1000000000000000000		     E-16 */
+      /*   lo:	  .................1000000000000000  E-16 */
       Int newexp=MINI(hi->exponent, hi->exponent+hilen-DECPMAX)-3;
-    Int reduce=newexp-lo->exponent;
-    if (reduce>0) {			/* [= case gives reduce=0 nop] */
+
       /* printf("FMA reduce: %ld\n", (LI)reduce); */
-      if (reduce>=lolen) {		/* eating all */
-	lo->lsd=lo->msd;		/* reduce to single digit */
-	lo->exponent=newexp;		/* [known to be non-zero] */
-	}
-       else { /* < */
-	uByte *up=lo->lsd;
-	lo->lsd=lo->lsd-reduce;
-	if (*lo->lsd==0)		/* could need sticky bit */
-	 for (; up>lo->lsd; up--) {	/* search discarded digits */
-	  if (*up!=0) {			/* found one... */
-	    *lo->lsd=1;			/* set sticky bit */
-	    break;
-	    }
-	  }
-	lo->exponent+=reduce;
-	}
+      lo->lsd=lo->msd;			     /* to single digit [maybe 0] */
+      lo->exponent=newexp;		     /* new lowest exponent */
       padding=hi->exponent-lo->exponent;     /* recalculate */
-      ul=lo->lsd;			 /* .. */
-      } /* maybe reduce */
-    /* padding is now <= DECPMAX+2 but still > 0; tricky DOUBLE case */
-    /* is when hi is a 1 that will become a 0.9999... by subtraction: */
-    /*	 hi:   1				   E+16 */
-    /*	 lo:	.................1000000000000000  E-16 */
-    /* which for the addition pads and reduces to: */
-    /*	 hi:   1000000000000000000		   E-2 */
-    /*	 lo:	.................1		   E-2 */
+      ul=lo->lsd;			     /* .. and repoint */
+      }
+
+    /* padding is still > 0, but will fit in acc (less leading carry slot) */
     #if DECCHECK
-      if (padding>DECPMAX+2) printf("FMA excess padding: %ld\n", (LI)padding);
       if (padding<=0) printf("FMA low padding: %ld\n", (LI)padding);
+      if (hilen+padding+1>FMALEN)
+	printf("FMA excess hilen+padding: %ld+%ld \n", (LI)hilen, (LI)padding);
       /* printf("FMA padding: %ld\n", (LI)padding); */
     #endif
+
     /* padding digits can now be set in the result; one or more of */
     /* these will come from lo; others will be zeros in the gap */
+    for (; ul-3>=lo->msd && padding>3; padding-=4, ul-=4, ub-=4) {
+      UBFROMUI(ub-3, UBTOUI(ul-3));	     /* [cannot overlap] */
+      }
     for (; ul>=lo->msd && padding>0; padding--, ul--, ub--) *ub=*ul;
     for (;padding>0; padding--, ub--) *ub=0; /* mind the gap */
     }
@@ -2088,23 +2196,39 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
   /* addition now complete to the right of the rightmost digit of hi */
   uh=hi->lsd;
 
-  /* carry was set up depending on ten's complement above; do the add... */
+  /* dow do the add from hi->lsd to the left */
+  /* [bytewise, because either operand can run out at any time] */
+  /* carry was set up depending on ten's complement above */
+  /* first assume both operands have some digits */
   for (;; ub--) {
-    uInt hid, lod;
-    if (uh<hi->msd) {
-      if (ul<lo->msd) break;
-      hid=hipad;
+    if (uh<hi->msd || ul<lo->msd) break;
+    *ub=(uByte)(carry+(*uh--)+(*ul--));
+    carry=0;
+    if (*ub<10) continue;
+    *ub-=10;
+    carry=1;
+    } /* both loop */
+
+  if (ul<lo->msd) {	      /* to left of lo */
+    for (;; ub--) {
+      if (uh<hi->msd) break;
+      *ub=(uByte)(carry+(*uh--));  /* [+0] */
+      carry=0;
+      if (*ub<10) continue;
+      *ub-=10;
+      carry=1;
+      } /* hi loop */
     }
-     else hid=*uh--;
-    if (ul<lo->msd) lod=0;
-     else lod=*ul--;
-    *ub=(uByte)(carry+hid+lod);
-    if (*ub<10) carry=0;
-     else {
+   else {		      /* to left of hi */
+    for (;; ub--) {
+      if (ul<lo->msd) break;
+      *ub=(uByte)(carry+hipad+(*ul--));
+      carry=0;
+      if (*ub<10) continue;
       *ub-=10;
       carry=1;
+      } /* lo loop */
     }
-    } /* addition loop */
 
   /* addition complete -- now handle carry, borrow, etc. */
   /* use lo to set up the num (its exponent is already correct, and */
@@ -2122,7 +2246,7 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
     if (!carry) {		   /* no carry out means hi<lo */
       /* borrowed -- take ten's complement of the right digits */
       lo->sign=hi->sign;	   /* sign is lhs sign */
-      for (ul=lo->msd; ul<lo->lsd-3; ul+=4) UINTAT(ul)=0x09090909-UINTAT(ul);
+      for (ul=lo->msd; ul<lo->lsd-3; ul+=4) UBFROMUI(ul, 0x09090909-UBTOUI(ul));
       for (; ul<=lo->lsd; ul++) *ul=(uByte)(0x09-*ul); /* [leaves ul at lsd+1] */
       /* complete the ten's complement by adding 1 [cannot overrun] */
       for (ul--; *ul==9; ul--) *ul=0;
@@ -2133,7 +2257,7 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
       /* all done except for the special IEEE 754 exact-zero-result */
       /* rule (see above); while testing for zero, strip leading */
       /* zeros (which will save decFinalize doing it) */
-      for (; UINTAT(lo->msd)==0 && lo->msd+3<lo->lsd;) lo->msd+=4;
+      for (; UBTOUI(lo->msd)==0 && lo->msd+3<lo->lsd;) lo->msd+=4;
       for (; *lo->msd==0 && lo->msd<lo->lsd;) lo->msd++;
       if (*lo->msd==0) {	   /* must be true zero (and diffsign) */
 	lo->sign=0;		   /* assume + */
@@ -2143,6 +2267,13 @@ decFloat * decFloatFMA(decFloat *result, const decFloat *dfl,
       } /* subtraction gave positive result */
     } /* diffsign */
 
+  #if DECCHECK
+  /* assert no left underrun */
+  if (lo->msd<acc) {
+    printf("FMA underrun by %ld \n", (LI)(acc-lo->msd));
+    }
+  #endif
+
   return decFinalize(result, lo, set);	/* round, check, and lay out */
   } /* decFloatFMA */
 
@@ -2242,11 +2373,11 @@ decFloat * decFloatInvert(decFloat *result, const decFloat *df,
 /* decFloatIs -- decFloat tests (IsSigned, etc.)		      */
 /*								      */
 /*   df is the decFloat to test 				      */
-/*   returns 0 or 1 in an int32_t				      */
+/*   returns 0 or 1 in a uInt					      */
 /*								      */
 /* Many of these could be macros, but having them as real functions   */
-/* is a bit cleaner (and they can be referred to here by the generic  */
-/* names)							      */
+/* is a little cleaner (and they can be referred to here by the       */
+/* generic names)						      */
 /* ------------------------------------------------------------------ */
 uInt decFloatIsCanonical(const decFloat *df) {
   if (DFISSPECIAL(df)) {
@@ -2333,7 +2464,7 @@ uInt decFloatIsZero(const decFloat *df) {
   } /* decFloatIs... */
 
 /* ------------------------------------------------------------------ */
-/* decFloatLogB -- return adjusted exponent, by 754r rules	      */
+/* decFloatLogB -- return adjusted exponent, by 754 rules	      */
 /*								      */
 /*   result gets the adjusted exponent as an integer, or a NaN etc.   */
 /*   df     is the decFloat to be examined			      */
@@ -2356,7 +2487,7 @@ decFloat * decFloatLogB(decFloat *result, const decFloat *df,
     return decInfinity(result, result);      /* canonical +Infinity */
     }
   if (DFISZERO(df)) {
-    set->status|=DEC_Division_by_zero;	     /* as per 754r */
+    set->status|=DEC_Division_by_zero;	     /* as per 754 */
     DFWORD(result, 0)=DECFLOAT_Sign;	     /* make negative */
     return decInfinity(result, result);      /* canonical -Infinity */
     }
@@ -2561,7 +2692,7 @@ decFloat * decFloatMultiply(decFloat *result,
 /*   set    is the context					      */
 /*   returns result						      */
 /*								      */
-/* This is 754r nextdown; Invalid is the only status possible (from   */
+/* This is 754 nextdown; Invalid is the only status possible (from    */
 /* an sNaN).							      */
 /* ------------------------------------------------------------------ */
 decFloat * decFloatNextMinus(decFloat *result, const decFloat *dfl,
@@ -2604,7 +2735,7 @@ decFloat * decFloatNextMinus(decFloat *result, const decFloat *dfl,
 /*   set    is the context					      */
 /*   returns result						      */
 /*								      */
-/* This is 754r nextup; Invalid is the only status possible (from     */
+/* This is 754 nextup; Invalid is the only status possible (from      */
 /* an sNaN).							      */
 /* ------------------------------------------------------------------ */
 decFloat * decFloatNextPlus(decFloat *result, const decFloat *dfl,
@@ -2649,8 +2780,9 @@ decFloat * decFloatNextPlus(decFloat *result, const decFloat *dfl,
 /*   set    is the context					      */
 /*   returns result						      */
 /*								      */
-/* This is 754r nextafter; status may be set unless the result is a   */
-/* normal number.						      */
+/* This is 754-1985 nextafter, as modified during revision (dropped   */
+/* from 754-2008); status may be set unless the result is a normal    */
+/* number.							      */
 /* ------------------------------------------------------------------ */
 decFloat * decFloatNextToward(decFloat *result,
 			      const decFloat *dfl, const decFloat *dfr,
@@ -2775,14 +2907,17 @@ decFloat * decFloatQuantize(decFloat *result,
 			    decContext *set) {
   Int	explb, exprb;	      /* left and right biased exponents */
   uByte *ulsd;		      /* local LSD pointer */
-  uInt	*ui;		      /* work */
-  uByte *ub;		      /* .. */
+  uByte *ub, *uc;	      /* work */
   Int	drop;		      /* .. */
   uInt	dpd;		      /* .. */
   uInt	encode; 	      /* encoding accumulator */
   uInt	sourhil, sourhir;     /* top words from source decFloats */
+  uInt	uiwork; 	      /* for macros */
+  #if QUAD
+  uShort uswork;	      /* .. */
+  #endif
   /* the following buffer holds the coefficient for manipulation */
-  uByte buf[4+DECPMAX*3];     /* + space for zeros to left or right */
+  uByte buf[4+DECPMAX*3+2*QUAD];   /* + space for zeros to left or right */
   #if DECTRACE
   bcdnum num;			   /* for trace displays */
   #endif
@@ -2838,7 +2973,7 @@ decFloat * decFloatQuantize(decFloat *result,
     /* there is at least one zero needed to the left, in all but one */
     /* exceptional (all-nines) case, so place four zeros now; this is */
     /* needed almost always and makes rounding all-nines by fours safe */
-    UINTAT(BUFOFF-4)=0;
+    UBFROMUI(BUFOFF-4, 0);
 
     /* Three cases here: */
     /*	 1. new LSD is in coefficient (almost always) */
@@ -2932,7 +3067,7 @@ decFloat * decFloatQuantize(decFloat *result,
 	/* increment the coefficient; this could give 1000... (after */
 	/* the all nines case) */
 	ub=ulsd;
-	for (; UINTAT(ub-3)==0x09090909; ub-=4) UINTAT(ub-3)=0;
+	for (; UBTOUI(ub-3)==0x09090909; ub-=4) UBFROMUI(ub-3, 0);
 	/* now at most 3 digits left to non-9 (usually just the one) */
 	for (; *ub==9; ub--) *ub=0;
 	*ub+=1;
@@ -2945,8 +3080,8 @@ decFloat * decFloatQuantize(decFloat *result,
     /* available in the coefficent -- the first word to the left was */
     /* cleared earlier for safe carry; now add any more needed */
     if (drop>4) {
-      UINTAT(BUFOFF-8)=0;		     /* must be at least 5 */
-      for (ui=&UINTAT(BUFOFF-12); ui>&UINTAT(ulsd-DECPMAX-3); ui--) *ui=0;
+      UBFROMUI(BUFOFF-8, 0);		     /* must be at least 5 */
+      for (uc=BUFOFF-12; uc>ulsd-DECPMAX-3; uc-=4) UBFROMUI(uc, 0);
       }
     } /* need round (drop>0) */
 
@@ -2967,18 +3102,21 @@ decFloat * decFloatQuantize(decFloat *result,
       #else
       static const uInt dmask[]={0, 0xff000000, 0xffff0000, 0xffffff00};
       #endif
-      for (ui=&UINTAT(BUFOFF+DECPMAX);; ui++) {
-	*ui=0;
-	if (UINTAT(&UBYTEAT(ui)-DECPMAX)!=0) { /* could be bad */
+      /* note that here zeros to the right are added by fours, so in */
+      /* the Quad case this could write 36 zeros if the coefficient has */
+      /* fewer than three significant digits (hence the +2*QUAD for buf) */
+      for (uc=BUFOFF+DECPMAX;; uc+=4) {
+	UBFROMUI(uc, 0);
+	if (UBTOUI(uc-DECPMAX)!=0) {		  /* could be bad */
 	  /* if all four digits should be zero, definitely bad */
-	  if (ui<=&UINTAT(BUFOFF+DECPMAX+(-drop)-4))
+	  if (uc<=BUFOFF+DECPMAX+(-drop)-4)
 	    return decInvalid(result, set);
 	  /* must be a 1- to 3-digit sequence; check more carefully */
-	  if ((UINTAT(&UBYTEAT(ui)-DECPMAX)&dmask[(-drop)%4])!=0)
+	  if ((UBTOUI(uc-DECPMAX)&dmask[(-drop)%4])!=0)
 	    return decInvalid(result, set);
 	  break;    /* no need for loop end test */
 	  }
-	if (ui>=&UINTAT(BUFOFF+DECPMAX+(-drop)-4)) break; /* done */
+	if (uc>=BUFOFF+DECPMAX+(-drop)-4) break;  /* done */
 	}
       ulsd=BUFOFF+DECPMAX+(-drop)-1;
       } /* pad and check leading zeros */
@@ -3204,7 +3342,7 @@ uInt decFloatSameQuantum(const decFloat *dfl, const decFloat *dfr) {
   } /* decFloatSameQuantum */
 
 /* ------------------------------------------------------------------ */
-/* decFloatScaleB -- multiply by a power of 10, as per 754r	      */
+/* decFloatScaleB -- multiply by a power of 10, as per 754	      */
 /*								      */
 /*   result gets the result of the operation			      */
 /*   dfl    is the first decFloat (lhs) 			      */
@@ -3270,6 +3408,7 @@ decFloat * decFloatShift(decFloat *result,
   uByte  buf[DECPMAX*2];		/* coefficient + padding */
   uInt	 digits, savestat;		/* work */
   bcdnum num;				/* .. */
+  uInt	 uiwork;			/* for macros */
 
   if (DFISNAN(dfl)||DFISNAN(dfr)) return decNaNs(result, dfl, dfr, set);
   if (!DFISINT(dfr)) return decInvalid(result, set);
@@ -3299,8 +3438,8 @@ decFloat * decFloatShift(decFloat *result,
     num.lsd=buf+DECPMAX-shift-1;
     }
    else { /* shift left -- zero padding needed to right */
-    UINTAT(buf+DECPMAX)=0;		/* 8 will handle most cases */
-    UINTAT(buf+DECPMAX+4)=0;		/* .. */
+    UBFROMUI(buf+DECPMAX, 0);		/* 8 will handle most cases */
+    UBFROMUI(buf+DECPMAX+4, 0); 	/* .. */
     if (shift>8) memset(buf+DECPMAX+8, 0, 8+QUAD*18); /* all other cases */
     num.msd+=shift;
     num.lsd=num.msd+DECPMAX-1;
@@ -3489,6 +3628,7 @@ static Int decNumCompare(const decFloat *dfl, const decFloat *dfr, Flag tot) {
   Int	sigl, sigr;			/* LHS and RHS non-0 signums */
   Int	shift;				/* shift needed to align operands */
   uByte *ub, *uc;			/* work */
+  uInt	uiwork; 			/* for macros */
   /* buffers +2 if Quad (36 digits), need double plus 4 for safe padding */
   uByte bufl[DECPMAX*2+QUAD*2+4];	/* for LHS coefficient + padding */
   uByte bufr[DECPMAX*2+QUAD*2+4];	/* for RHS coefficient + padding */
@@ -3544,17 +3684,16 @@ static Int decNumCompare(const decFloat *dfl, const decFloat *dfr, Flag tot) {
   /* decode the coefficients */
   /* (shift both right two if Quad to make a multiple of four) */
   #if QUAD
-    ub=bufl;                            /* avoid type-pun violation */
-    UINTAT(ub)=0;
-    uc=bufr;                            /* avoid type-pun violation */
-    UINTAT(uc)=0;
+    UBFROMUI(bufl, 0);
+    UBFROMUI(bufr, 0);
   #endif
   GETCOEFF(dfl, bufl+QUAD*2);		/* decode from decFloat */
   GETCOEFF(dfr, bufr+QUAD*2);		/* .. */
   if (shift==0) {			/* aligned; common and easy */
     /* all multiples of four, here */
     for (ub=bufl, uc=bufr; ub<bufl+DECPMAX+QUAD*2; ub+=4, uc+=4) {
-      if (UINTAT(ub)==UINTAT(uc)) continue; /* so far so same */
+      uInt ui=UBTOUI(ub);
+      if (ui==UBTOUI(uc)) continue;	/* so far so same */
       /* about to find a winner; go by bytes in case little-endian */
       for (;; ub++, uc++) {
 	if (*ub>*uc) return sigl;	/* difference found */
@@ -3565,17 +3704,17 @@ static Int decNumCompare(const decFloat *dfl, const decFloat *dfr, Flag tot) {
    else if (shift>0) {			/* lhs to left */
     ub=bufl;				/* RHS pointer */
     /* pad bufl so right-aligned; most shifts will fit in 8 */
-    UINTAT(bufl+DECPMAX+QUAD*2)=0;	/* add eight zeros */
-    UINTAT(bufl+DECPMAX+QUAD*2+4)=0;	/* .. */
+    UBFROMUI(bufl+DECPMAX+QUAD*2, 0);	/* add eight zeros */
+    UBFROMUI(bufl+DECPMAX+QUAD*2+4, 0); /* .. */
     if (shift>8) {
       /* more than eight; fill the rest, and also worth doing the */
       /* lead-in by fours */
       uByte *up;			/* work */
       uByte *upend=bufl+DECPMAX+QUAD*2+shift;
-      for (up=bufl+DECPMAX+QUAD*2+8; up<upend; up+=4) UINTAT(up)=0;
+      for (up=bufl+DECPMAX+QUAD*2+8; up<upend; up+=4) UBFROMUI(up, 0);
       /* [pads up to 36 in all for Quad] */
       for (;; ub+=4) {
-	if (UINTAT(ub)!=0) return sigl;
+	if (UBTOUI(ub)!=0) return sigl;
 	if (ub+4>bufl+shift-4) break;
 	}
       }
@@ -3585,7 +3724,8 @@ static Int decNumCompare(const decFloat *dfl, const decFloat *dfr, Flag tot) {
     /* comparison can go for the full length of bufr, which is a */
     /* multiple of 4 bytes */
     for (uc=bufr; ; uc+=4, ub+=4) {
-      if (UINTAT(uc)!=UINTAT(ub)) {	/* mismatch found */
+      uInt ui=UBTOUI(ub);
+      if (ui!=UBTOUI(uc)) {		/* mismatch found */
 	for (;; uc++, ub++) {		/* check from left [little-endian?] */
 	  if (*ub>*uc) return sigl;	/* difference found */
 	  if (*ub<*uc) return sigr;	/* .. */
@@ -3598,17 +3738,17 @@ static Int decNumCompare(const decFloat *dfl, const decFloat *dfr, Flag tot) {
    else { /* shift<0) .. RHS is to left of LHS; mirror shift>0 */
     uc=bufr;				/* RHS pointer */
     /* pad bufr so right-aligned; most shifts will fit in 8 */
-    UINTAT(bufr+DECPMAX+QUAD*2)=0;	/* add eight zeros */
-    UINTAT(bufr+DECPMAX+QUAD*2+4)=0;	/* .. */
+    UBFROMUI(bufr+DECPMAX+QUAD*2, 0);	/* add eight zeros */
+    UBFROMUI(bufr+DECPMAX+QUAD*2+4, 0); /* .. */
     if (shift<-8) {
       /* more than eight; fill the rest, and also worth doing the */
       /* lead-in by fours */
       uByte *up;			/* work */
       uByte *upend=bufr+DECPMAX+QUAD*2-shift;
-      for (up=bufr+DECPMAX+QUAD*2+8; up<upend; up+=4) UINTAT(up)=0;
+      for (up=bufr+DECPMAX+QUAD*2+8; up<upend; up+=4) UBFROMUI(up, 0);
       /* [pads up to 36 in all for Quad] */
       for (;; uc+=4) {
-	if (UINTAT(uc)!=0) return sigr;
+	if (UBTOUI(uc)!=0) return sigr;
 	if (uc+4>bufr-shift-4) break;
 	}
       }
@@ -3618,7 +3758,8 @@ static Int decNumCompare(const decFloat *dfl, const decFloat *dfr, Flag tot) {
     /* comparison can go for the full length of bufl, which is a */
     /* multiple of 4 bytes */
     for (ub=bufl; ; ub+=4, uc+=4) {
-      if (UINTAT(ub)!=UINTAT(uc)) {	/* mismatch found */
+      uInt ui=UBTOUI(ub);
+      if (ui!=UBTOUI(uc)) {		/* mismatch found */
 	for (;; ub++, uc++) {		/* check from left [little-endian?] */
 	  if (*ub>*uc) return sigl;	/* difference found */
 	  if (*ub<*uc) return sigr;	/* .. */

libdecnumber/decCommon.c

@@ -161,7 +161,7 @@ static const uInt DECCOMBWEXP[64]={
 #if QUAD
 const uInt DECCOMBMSD[64]={
   0, 1, 2, 3, 4, 5, 6, 7, 0, 1, 2, 3, 4, 5, 6, 7,
-  0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 9, 8, 9, 0, 1,
+  0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 9, 8, 9, 0, 0,
   0, 1, 2, 3, 4, 5, 6, 7, 0, 1, 2, 3, 4, 5, 6, 7,
   0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 9, 8, 9, 0, 0};
 
@@ -248,6 +248,7 @@ static decFloat * decFinalize(decFloat *df, bcdnum *num,
 			      decContext *set) {
   uByte *ub;		      /* work */
   uInt	 dpd;		      /* .. */
+  uInt	 uiwork;	      /* for macros */
   uByte *umsd=num->msd;       /* local copy */
   uByte *ulsd=num->lsd;       /* .. */
   uInt	 encode;	      /* encoding accumulator */
@@ -279,7 +280,7 @@ static decFloat * decFinalize(decFloat *df, bcdnum *num,
     /* skip leading insignificant zeros to calculate an exact length */
     /* [this is quite expensive] */
     if (*umsd==0) {
-      for (; UINTAT(umsd)==0 && umsd+3<ulsd;) umsd+=4;
+      for (; umsd+3<ulsd && UBTOUI(umsd)==0;) umsd+=4;
       for (; *umsd==0 && umsd<ulsd;) umsd++;
       length=ulsd-umsd+1;		     /* recalculate */
       }
@@ -388,7 +389,9 @@ static decFloat * decFinalize(decFloat *df, bcdnum *num,
 	  /* increment the coefficient; this might end up with 1000... */
 	  /* (after the all nines case) */
 	  ub=ulsd;
-	  for(; ub-3>=umsd && UINTAT(ub-3)==0x09090909; ub-=4) UINTAT(ub-3)=0;
+	  for(; ub-3>=umsd && UBTOUI(ub-3)==0x09090909; ub-=4)	{
+	    UBFROMUI(ub-3, 0);		     /* to 00000000 */
+	    }
 	  /* [note ub could now be to left of msd, and it is not safe */
 	  /* to write to the the left of the msd] */
 	  /* now at most 3 digits left to non-9 (usually just the one) */
@@ -475,8 +478,8 @@ static decFloat * decFinalize(decFloat *df, bcdnum *num,
 	  uByte *t=buffer;		/* safe target */
 	  uByte *tlsd=buffer+(ulsd-umsd)+shift; /* target LSD */
 	  /* printf("folddown shift=%ld\n", (LI)shift); */
-	  for (; s<=ulsd; s+=4, t+=4) UINTAT(t)=UINTAT(s);
-	  for (t=tlsd-shift+1; t<=tlsd; t+=4) UINTAT(t)=0;  /* pad */
+	  for (; s<=ulsd; s+=4, t+=4) UBFROMUI(t, UBTOUI(s));
+	  for (t=tlsd-shift+1; t<=tlsd; t+=4) UBFROMUI(t, 0);  /* pad 0s */
 	  num->exponent-=shift;
 	  umsd=buffer;
 	  ulsd=tlsd;
@@ -492,23 +495,23 @@ static decFloat * decFinalize(decFloat *df, bcdnum *num,
   /*------------------------------------------------------------------*/
   /* Following code does not alter coefficient (could be allnines array) */
 
+  /* fast path possible when DECPMAX digits */
   if (length==DECPMAX) {
     return decFloatFromBCD(df, num->exponent, umsd, num->sign);
-    }
+    } /* full-length */
 
-  /* Here when length is short */
+  /* slower path when not a full-length number; must care about length */
+  /* [coefficient length here will be < DECPMAX] */
   if (!NUMISSPECIAL(num)) {		/* is still finite */
     /* encode the combination field and exponent continuation */
     uInt uexp=(uInt)(num->exponent+DECBIAS); /* biased exponent */
     uInt code=(uexp>>DECECONL)<<4;	/* top two bits of exp */
-    /* [msd=0] */
+    /* [msd==0] */
     /* look up the combination field and make high word */
     encode=DECCOMBFROM[code];		/* indexed by (0-2)*16+msd */
     encode|=(uexp<<(32-6-DECECONL)) & 0x03ffffff; /* exponent continuation */
     }
    else encode=num->exponent;		/* special [already in word] */
-  /* [coefficient length here will be < DECPMAX] */
-
   encode|=num->sign;			/* add sign */
 
   /* private macro to extract a declet, n (where 0<=n<DECLETS and 0 */
@@ -519,7 +522,7 @@ static decFloat * decFinalize(decFloat *df, bcdnum *num,
   /* working pointer, uInt *ub. */
   /* As not full-length then chances are there are many leading zeros */
   /* [and there may be a partial triad] */
-  #define getDPD(dpd, n) ub=ulsd-(3*(n))-2;			      \
+  #define getDPDt(dpd, n) ub=ulsd-(3*(n))-2;			      \
     if (ub<umsd-2) dpd=0;					      \
      else if (ub>=umsd) dpd=BCD2DPD[(*ub*256)+(*(ub+1)*16)+*(ub+2)];  \
      else {dpd=*(ub+2); if (ub+1==umsd) dpd+=*(ub+1)*16; dpd=BCD2DPD[dpd];}
@@ -528,48 +531,48 @@ static decFloat * decFinalize(decFloat *df, bcdnum *num,
   /* according to endianness; in all cases complete the sign word */
   /* first */
   #if DECPMAX==7
-    getDPD(dpd, 1);
+    getDPDt(dpd, 1);
     encode|=dpd<<10;
-    getDPD(dpd, 0);
+    getDPDt(dpd, 0);
     encode|=dpd;
     DFWORD(df, 0)=encode;     /* just the one word */
 
   #elif DECPMAX==16
-    getDPD(dpd, 4); encode|=dpd<<8;
-    getDPD(dpd, 3); encode|=dpd>>2;
+    getDPDt(dpd, 4); encode|=dpd<<8;
+    getDPDt(dpd, 3); encode|=dpd>>2;
     DFWORD(df, 0)=encode;
     encode=dpd<<30;
-    getDPD(dpd, 2); encode|=dpd<<20;
-    getDPD(dpd, 1); encode|=dpd<<10;
-    getDPD(dpd, 0); encode|=dpd;
+    getDPDt(dpd, 2); encode|=dpd<<20;
+    getDPDt(dpd, 1); encode|=dpd<<10;
+    getDPDt(dpd, 0); encode|=dpd;
     DFWORD(df, 1)=encode;
 
   #elif DECPMAX==34
-    getDPD(dpd,10); encode|=dpd<<4;
-    getDPD(dpd, 9); encode|=dpd>>6;
+    getDPDt(dpd,10); encode|=dpd<<4;
+    getDPDt(dpd, 9); encode|=dpd>>6;
     DFWORD(df, 0)=encode;
 
     encode=dpd<<26;
-    getDPD(dpd, 8); encode|=dpd<<16;
-    getDPD(dpd, 7); encode|=dpd<<6;
-    getDPD(dpd, 6); encode|=dpd>>4;
+    getDPDt(dpd, 8); encode|=dpd<<16;
+    getDPDt(dpd, 7); encode|=dpd<<6;
+    getDPDt(dpd, 6); encode|=dpd>>4;
     DFWORD(df, 1)=encode;
 
     encode=dpd<<28;
-    getDPD(dpd, 5); encode|=dpd<<18;
-    getDPD(dpd, 4); encode|=dpd<<8;
-    getDPD(dpd, 3); encode|=dpd>>2;
+    getDPDt(dpd, 5); encode|=dpd<<18;
+    getDPDt(dpd, 4); encode|=dpd<<8;
+    getDPDt(dpd, 3); encode|=dpd>>2;
     DFWORD(df, 2)=encode;
 
     encode=dpd<<30;
-    getDPD(dpd, 2); encode|=dpd<<20;
-    getDPD(dpd, 1); encode|=dpd<<10;
-    getDPD(dpd, 0); encode|=dpd;
+    getDPDt(dpd, 2); encode|=dpd<<20;
+    getDPDt(dpd, 1); encode|=dpd<<10;
+    getDPDt(dpd, 0); encode|=dpd;
     DFWORD(df, 3)=encode;
   #endif
 
   /* printf("Status: %08lx\n", (LI)set->status); */
-  /* decFloatShow(df, "final"); */
+  /* decFloatShow(df, "final2"); */
   return df;
   } /* decFinalize */
 
@@ -582,7 +585,7 @@ static decFloat * decFinalize(decFloat *df, bcdnum *num,
 /*  bcdar holds DECPMAX digits to set the coefficient from, one       */
 /*    digit in each byte (BCD8 encoding); the first (MSD) is ignored  */
 /*    if df is a NaN; all are ignored if df is infinite.	      */
-/*    All bytes must be in 0-9; results undefined otherwise.	      */
+/*    All bytes must be in 0-9; results are undefined otherwise.      */
 /*  sig is DECFLOAT_Sign to set the sign bit, 0 otherwise	      */
 /*  returns df, which will be canonical 			      */
 /*								      */
@@ -609,53 +612,53 @@ decFloat * decFloatFromBCD(decFloat *df, Int exp, const uByte *bcdar,
   /* and put the corresponding DPD code into dpd. */
   /* Use of a working pointer, uInt *ub, is assumed. */
 
-  #define getDPDf(dpd, n) ub=bcdar+DECPMAX-1-(3*(n))-2;	    \
+  #define getDPDb(dpd, n) ub=bcdar+DECPMAX-1-(3*(n))-2;     \
     dpd=BCD2DPD[(*ub*256)+(*(ub+1)*16)+*(ub+2)];
 
   /* place the declets in the encoding words and copy to result (df), */
   /* according to endianness; in all cases complete the sign word */
   /* first */
   #if DECPMAX==7
-    getDPDf(dpd, 1);
+    getDPDb(dpd, 1);
     encode|=dpd<<10;
-    getDPDf(dpd, 0);
+    getDPDb(dpd, 0);
     encode|=dpd;
     DFWORD(df, 0)=encode;     /* just the one word */
 
   #elif DECPMAX==16
-    getDPDf(dpd, 4); encode|=dpd<<8;
-    getDPDf(dpd, 3); encode|=dpd>>2;
+    getDPDb(dpd, 4); encode|=dpd<<8;
+    getDPDb(dpd, 3); encode|=dpd>>2;
     DFWORD(df, 0)=encode;
     encode=dpd<<30;
-    getDPDf(dpd, 2); encode|=dpd<<20;
-    getDPDf(dpd, 1); encode|=dpd<<10;
-    getDPDf(dpd, 0); encode|=dpd;
+    getDPDb(dpd, 2); encode|=dpd<<20;
+    getDPDb(dpd, 1); encode|=dpd<<10;
+    getDPDb(dpd, 0); encode|=dpd;
     DFWORD(df, 1)=encode;
 
   #elif DECPMAX==34
-    getDPDf(dpd,10); encode|=dpd<<4;
-    getDPDf(dpd, 9); encode|=dpd>>6;
+    getDPDb(dpd,10); encode|=dpd<<4;
+    getDPDb(dpd, 9); encode|=dpd>>6;
     DFWORD(df, 0)=encode;
 
     encode=dpd<<26;
-    getDPDf(dpd, 8); encode|=dpd<<16;
-    getDPDf(dpd, 7); encode|=dpd<<6;
-    getDPDf(dpd, 6); encode|=dpd>>4;
+    getDPDb(dpd, 8); encode|=dpd<<16;
+    getDPDb(dpd, 7); encode|=dpd<<6;
+    getDPDb(dpd, 6); encode|=dpd>>4;
     DFWORD(df, 1)=encode;
 
     encode=dpd<<28;
-    getDPDf(dpd, 5); encode|=dpd<<18;
-    getDPDf(dpd, 4); encode|=dpd<<8;
-    getDPDf(dpd, 3); encode|=dpd>>2;
+    getDPDb(dpd, 5); encode|=dpd<<18;
+    getDPDb(dpd, 4); encode|=dpd<<8;
+    getDPDb(dpd, 3); encode|=dpd>>2;
     DFWORD(df, 2)=encode;
 
     encode=dpd<<30;
-    getDPDf(dpd, 2); encode|=dpd<<20;
-    getDPDf(dpd, 1); encode|=dpd<<10;
-    getDPDf(dpd, 0); encode|=dpd;
+    getDPDb(dpd, 2); encode|=dpd<<20;
+    getDPDb(dpd, 1); encode|=dpd<<10;
+    getDPDb(dpd, 0); encode|=dpd;
     DFWORD(df, 3)=encode;
   #endif
-  /* decFloatShow(df, "final"); */
+  /* decFloatShow(df, "fromB"); */
   return df;
   } /* decFloatFromBCD */
 
@@ -702,6 +705,70 @@ decFloat * decFloatFromPacked(decFloat *df, Int exp, const uByte *packed) {
   } /* decFloatFromPacked */
 
 /* ------------------------------------------------------------------ */
+/* decFloatFromPackedChecked -- set from exponent and packed; checked */
+/*								      */
+/*  df is the target decFloat					      */
+/*  exp is the in-range unbiased exponent, q, or a special value in   */
+/*    the form returned by decFloatGetExponent			      */
+/*  packed holds DECPMAX packed decimal digits plus a sign nibble     */
+/*    (all 6 codes are OK); the first (MSD) must be 0 if df is a NaN  */
+/*    and all digits must be 0 if df is infinite.  For DOUBLE and     */
+/*    QUAD the first (pad) nibble must be 0.			      */
+/*    All coefficient nibbles must be in 0-9 and sign in A-F.	      */
+/*  returns df, which will be canonical or NULL if any of the	      */
+/*    requirements are not met (if this case df is unchanged); that   */
+/*    is, the input data must be as returned by decFloatToPacked,     */
+/*    except that all six sign codes are acccepted.		      */
+/*								      */
+/* No status will be set.					      */
+/* ------------------------------------------------------------------ */
+decFloat * decFloatFromPackedChecked(decFloat *df, Int exp,
+				     const uByte *packed) {
+  uByte bcdar[DECPMAX+2];		/* work [+1 for pad, +1 for sign] */
+  const uByte *ip;			/* .. */
+  uByte *op;				/* .. */
+  Int	sig=0;				/* sign */
+
+  /* expand coefficient and sign to BCDAR */
+  #if SINGLE
+  op=bcdar+1;				/* no pad digit */
+  #else
+  op=bcdar;				/* first (pad) digit here */
+  #endif
+  for (ip=packed; ip<packed+((DECPMAX+2)/2); ip++) {
+    *op=*ip>>4;
+    if (*op>9) return NULL;
+    op++;
+    *op=(uByte)(*ip&0x0f);		/* [final nibble is sign] */
+    if (*op>9 && ip<packed+((DECPMAX+2)/2)-1) return NULL;
+    op++;
+    }
+  op--; 				/* -> sign byte */
+  if (*op<=9) return NULL;		/* bad sign */
+  if (*op==DECPMINUS || *op==DECPMINUSALT) sig=DECFLOAT_Sign;
+
+  #if !SINGLE
+  if (bcdar[0]!=0) return NULL; 	/* bad pad nibble */
+  #endif
+
+  if (EXPISNAN(exp)) {			/* a NaN */
+    if (bcdar[1]!=0) return NULL;	/* bad msd */
+    } /* NaN */
+   else if (EXPISINF(exp)) {		/* is infinite */
+    Int i;
+    for (i=0; i<DECPMAX; i++) {
+      if (bcdar[i+1]!=0) return NULL;	/* should be all zeros */
+      }
+    } /* infinity */
+   else {				/* finite */
+    /* check the exponent is in range */
+    if (exp>DECEMAX-DECPMAX+1) return NULL;
+    if (exp<DECEMIN-DECPMAX+1) return NULL;
+    }
+  return decFloatFromBCD(df, exp, bcdar+1, sig);
+  } /* decFloatFromPacked */
+
+/* ------------------------------------------------------------------ */
 /* decFloatFromString -- conversion from numeric string 	      */
 /*								      */
 /*  result  is the decFloat format number which gets the result of    */
@@ -730,6 +797,7 @@ decFloat * decFloatFromString(decFloat *result, const char *string,
   const  char *cfirst=string;	   /* -> first character of decimal part */
   const  char *c;		   /* work */
   uByte *ub;			   /* .. */
+  uInt	 uiwork;		   /* for macros */
   bcdnum num;			   /* collects data for finishing */
   uInt	 error=DEC_Conversion_syntax;	/* assume the worst */
   uByte  buffer[ROUNDUP(DECSTRING+11, 8)]; /* room for most coefficents, */
@@ -832,8 +900,8 @@ decFloat * decFloatFromString(decFloat *result, const char *string,
 	    /* as usual, go by fours when safe; NB it has been asserted */
 	    /* that a '.' does not have the same mask as a digit */
 	    if (c<=clast-3			       /* safe for four */
-	     && (UINTAT(c)&0xf0f0f0f0)==CHARMASK) {    /* test four */
-	      UINTAT(ub)=UINTAT(c)&0x0f0f0f0f;	       /* to BCD8 */
+	     && (UBTOUI(c)&0xf0f0f0f0)==CHARMASK) {    /* test four */
+	      UBFROMUI(ub, UBTOUI(c)&0x0f0f0f0f);      /* to BCD8 */
 	      ub+=4;
 	      c+=4;
 	      continue;
@@ -846,7 +914,7 @@ decFloat * decFloatFromString(decFloat *result, const char *string,
 	    }
 	  } /* had dot */
 	/* Now no dot; do this by fours (where safe) */
-	for (; c<=clast-3; c+=4, ub+=4) UINTAT(ub)=UINTAT(c)&0x0f0f0f0f;
+	for (; c<=clast-3; c+=4, ub+=4) UBFROMUI(ub, UBTOUI(c)&0x0f0f0f0f);
 	for (; c<=clast; c++, ub++) *ub=(uByte)(*c-'0');
 	num.lsd=buffer+digits-1;	     /* record new LSD */
 	} /* fits */
@@ -871,8 +939,8 @@ decFloat * decFloatFromString(decFloat *result, const char *string,
 	for (c=cfirst; c<=clast && ub<=buffer+DECPMAX; c++) {
 	  /* (see commentary just above) */
 	  if (c<=clast-3			  /* safe for four */
-	   && (UINTAT(c)&0xf0f0f0f0)==CHARMASK) { /* four digits */
-	    UINTAT(ub)=UINTAT(c)&0x0f0f0f0f;	  /* to BCD8 */
+	   && (UBTOUI(c)&0xf0f0f0f0)==CHARMASK) { /* four digits */
+	    UBFROMUI(ub, UBTOUI(c)&0x0f0f0f0f);   /* to BCD8 */
 	    ub+=4;
 	    c+=3;			     /* [will become 4] */
 	    continue;
@@ -1174,6 +1242,9 @@ char * decFloatToEngString(const decFloat *df, char *string){
   char *s, *t;			   /* .. (source, target) */
   Int  pre, e;			   /* work */
   const uByte *u;		   /* .. */
+  uInt	uiwork; 		   /* for macros [one compiler needs */
+				   /* volatile here to avoid bug, but */
+				   /* that doubles execution time] */
 
   /* Source words; macro handles endianness */
   uInt sourhi=DFWORD(df, 0);	   /* word with sign */
@@ -1231,9 +1302,9 @@ char * decFloatToEngString(const decFloat *df, char *string){
   /* and are safe because the last item in the array is of length */
   /* three and has the length byte following.) */
   #define dpd2char(dpdin) u=&DPD2BCD8[((dpdin)&0x3ff)*4];	 \
-	 if (c!=cstart) {UINTAT(c)=UINTAT(u)|CHARMASK; c+=3;}	 \
+	 if (c!=cstart) {UBFROMUI(c, UBTOUI(u)|CHARMASK); c+=3;} \
 	  else if (*(u+3)) {					 \
-	   UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; c+=*(u+3);}
+	   UBFROMUI(c, UBTOUI(u+3-*(u+3))|CHARMASK); c+=*(u+3);}
 
   #if DECPMAX==7
   dpd2char(sourhi>>10); 		/* declet 1 */
@@ -1310,8 +1381,8 @@ char * decFloatToEngString(const decFloat *df, char *string){
       /* because there is still space for exponent */
       s=dotat+ROUNDDOWN4(c-dotat);	/* source */
       t=s+1;				/* target */
-      /* open the gap */
-      for (; s>=dotat; s-=4, t-=4) UINTAT(t)=UINTAT(s);
+      /* open the gap [cannot use memcpy] */
+      for (; s>=dotat; s-=4, t-=4) UBFROMUI(t, UBTOUI(s));
       *dotat='.';
       c++;				/* length increased by one */
       } /* need dot? */
@@ -1326,19 +1397,19 @@ char * decFloatToEngString(const decFloat *df, char *string){
     /* backoff if too far to the right */
     if (t>string+DECSTRING-5) t=string+DECSTRING-5; /* adjust to fit */
     /* now shift the entire coefficient to the right, being careful not */
-    /* to access to the left of string */
-    for (s=t-pre; s>=string; s-=4, t-=4) UINTAT(t)=UINTAT(s);
+    /* to access to the left of string [cannot use memcpy] */
+    for (s=t-pre; s>=string; s-=4, t-=4) UBFROMUI(t, UBTOUI(s));
     /* for Quads and Singles there may be a character or two left... */
     s+=3;				/* where next would come from */
     for(; s>=cstart; s--, t--) *(t+3)=*(s);
     /* now have fill 0. through 0.00000; use overlaps to avoid tests */
     if (pre>=4) {
-      UINTAT(cstart+pre-4)=UINTAT("0000");
-      UINTAT(cstart)=UINTAT("0.00");
+      memcpy(cstart+pre-4, "0000", 4);
+      memcpy(cstart, "0.00", 4);
       }
      else { /* 2 or 3 */
       *(cstart+pre-1)='0';
-      USHORTAT(cstart)=USHORTAT("0.");
+      memcpy(cstart, "0.", 2);
       }
     c+=pre;				/* to end */
     }
@@ -1346,7 +1417,7 @@ char * decFloatToEngString(const decFloat *df, char *string){
   /* finally add the E-part, if needed; it will never be 0, and has */
   /* a maximum length of 3 or 4 digits (asserted above) */
   if (e!=0) {
-    USHORTAT(c)=USHORTAT("E+");		/* starts with E, assume + */
+    memcpy(c, "E+", 2); 		/* starts with E, assume + */
     c++;
     if (e<0) {
       *c='-';				/* oops, need '-' */
@@ -1358,12 +1429,12 @@ char * decFloatToEngString(const decFloat *df, char *string){
       u=&BIN2BCD8[e*4]; 		/* -> 3 digits + length byte */
       /* copy fixed 4 characters [is safe], starting at non-zero */
       /* and with character mask to convert BCD to char */
-      UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK;
+      UBFROMUI(c, UBTOUI(u+3-*(u+3))|CHARMASK);
       c+=*(u+3);			/* bump pointer appropriately */
     #elif DECEMAXD==4
       if (e<1000) {			/* 3 (or fewer) digits case */
 	u=&BIN2BCD8[e*4];		/* -> 3 digits + length byte */
-	UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; /* [as above] */
+	UBFROMUI(c, UBTOUI(u+3-*(u+3))|CHARMASK); /* [as above] */
 	c+=*(u+3);			/* bump pointer appropriately */
 	}
        else {				/* 4-digits */
@@ -1371,7 +1442,7 @@ char * decFloatToEngString(const decFloat *df, char *string){
 	Int rem=e-(1000*thou);		/* e%1000 */
 	*c++=(char)('0'+(char)thou);	/* the thousands digit */
 	u=&BIN2BCD8[rem*4];		/* -> 3 digits + length byte */
-	UINTAT(c)=UINTAT(u)|CHARMASK;	/* copy fixed 3+1 characters [is safe] */
+	UBFROMUI(c, UBTOUI(u)|CHARMASK);/* copy fixed 3+1 characters [is safe] */
 	c+=3;				/* bump pointer, always 3 digits */
 	}
     #endif
@@ -1450,6 +1521,9 @@ char * decFloatToString(const decFloat *df, char *string){
   char *s, *t;			   /* .. (source, target) */
   Int  pre, e;			   /* work */
   const uByte *u;		   /* .. */
+  uInt	uiwork; 		   /* for macros [one compiler needs */
+				   /* volatile here to avoid bug, but */
+				   /* that doubles execution time] */
 
   /* Source words; macro handles endianness */
   uInt sourhi=DFWORD(df, 0);	   /* word with sign */
@@ -1467,7 +1541,11 @@ char * decFloatToString(const decFloat *df, char *string){
   msd=DECCOMBMSD[comb]; 	   /* decode the combination field */
   exp=DECCOMBEXP[comb]; 	   /* .. */
 
-  if (EXPISSPECIAL(exp)) {	   /* special */
+  if (!EXPISSPECIAL(exp)) {	   /* finite */
+    /* complete exponent; top two bits are in place */
+    exp+=GETECON(df)-DECBIAS;	   /* .. + continuation and unbias */
+    }
+   else {			   /* IS special */
     if (exp==DECFLOAT_Inf) {	   /* infinity */
       strcpy(c, "Infinity");
       return string;		   /* easy */
@@ -1487,9 +1565,6 @@ char * decFloatToString(const decFloat *df, char *string){
     /* otherwise drop through to add integer; set correct exp etc. */
     exp=0; msd=0;		   /* setup for following code */
     }
-   else { /* complete exponent; top two bits are in place */
-    exp+=GETECON(df)-DECBIAS;	   /* .. + continuation and unbias */
-    }
 
   /* convert the digits of the significand to characters */
   cstart=c;			   /* save start of coefficient */
@@ -1506,9 +1581,9 @@ char * decFloatToString(const decFloat *df, char *string){
   /* and are safe because the last item in the array is of length */
   /* three and has the length byte following.) */
   #define dpd2char(dpdin) u=&DPD2BCD8[((dpdin)&0x3ff)*4];	 \
-	 if (c!=cstart) {UINTAT(c)=UINTAT(u)|CHARMASK; c+=3;}	 \
+	 if (c!=cstart) {UBFROMUI(c, UBTOUI(u)|CHARMASK); c+=3;} \
 	  else if (*(u+3)) {					 \
-	   UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; c+=*(u+3);}
+	   UBFROMUI(c, UBTOUI(u+3-*(u+3))|CHARMASK); c+=*(u+3);}
 
   #if DECPMAX==7
   dpd2char(sourhi>>10); 		/* declet 1 */
@@ -1556,12 +1631,13 @@ char * decFloatToString(const decFloat *df, char *string){
   if (pre>0) {			   /* ddd.ddd (plain), perhaps with E */
     char *dotat=cstart+pre;
     if (dotat<c) {			/* if embedded dot needed... */
+      /* [memmove is a disaster, here] */
       /* move by fours; there must be space for junk at the end */
-      /* because there is still space for exponent */
+      /* because exponent is still possible */
       s=dotat+ROUNDDOWN4(c-dotat);	/* source */
       t=s+1;				/* target */
-      /* open the gap */
-      for (; s>=dotat; s-=4, t-=4) UINTAT(t)=UINTAT(s);
+      /* open the gap [cannot use memcpy] */
+      for (; s>=dotat; s-=4, t-=4) UBFROMUI(t, UBTOUI(s));
       *dotat='.';
       c++;				/* length increased by one */
       } /* need dot? */
@@ -1569,7 +1645,7 @@ char * decFloatToString(const decFloat *df, char *string){
     /* finally add the E-part, if needed; it will never be 0, and has */
     /* a maximum length of 3 or 4 digits (asserted above) */
     if (e!=0) {
-      USHORTAT(c)=USHORTAT("E+");	/* starts with E, assume + */
+      memcpy(c, "E+", 2);		/* starts with E, assume + */
       c++;
       if (e<0) {
 	*c='-'; 			/* oops, need '-' */
@@ -1581,12 +1657,12 @@ char * decFloatToString(const decFloat *df, char *string){
 	u=&BIN2BCD8[e*4];		/* -> 3 digits + length byte */
 	/* copy fixed 4 characters [is safe], starting at non-zero */
 	/* and with character mask to convert BCD to char */
-	UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK;
+	UBFROMUI(c, UBTOUI(u+3-*(u+3))|CHARMASK);
 	c+=*(u+3);			/* bump pointer appropriately */
       #elif DECEMAXD==4
 	if (e<1000) {			/* 3 (or fewer) digits case */
 	  u=&BIN2BCD8[e*4];		/* -> 3 digits + length byte */
-	  UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; /* [as above] */
+	  UBFROMUI(c, UBTOUI(u+3-*(u+3))|CHARMASK); /* [as above] */
 	  c+=*(u+3);			/* bump pointer appropriately */
 	  }
 	 else { 			/* 4-digits */
@@ -1594,7 +1670,7 @@ char * decFloatToString(const decFloat *df, char *string){
 	  Int rem=e-(1000*thou);	/* e%1000 */
 	  *c++=(char)('0'+(char)thou);	/* the thousands digit */
 	  u=&BIN2BCD8[rem*4];		/* -> 3 digits + length byte */
-	  UINTAT(c)=UINTAT(u)|CHARMASK; /* copy fixed 3+1 characters [is safe] */
+	  UBFROMUI(c, UBTOUI(u)|CHARMASK); /* copy fixed 3+1 characters [is safe] */
 	  c+=3; 			/* bump pointer, always 3 digits */
 	  }
       #endif
@@ -1618,19 +1694,19 @@ char * decFloatToString(const decFloat *df, char *string){
   /* backoff if too far to the right */
   if (t>string+DECSTRING-5) t=string+DECSTRING-5; /* adjust to fit */
   /* now shift the entire coefficient to the right, being careful not */
-  /* to access to the left of string */
-  for (s=t-pre; s>=string; s-=4, t-=4) UINTAT(t)=UINTAT(s);
+  /* to access to the left of string [cannot use memcpy] */
+  for (s=t-pre; s>=string; s-=4, t-=4) UBFROMUI(t, UBTOUI(s));
   /* for Quads and Singles there may be a character or two left... */
   s+=3; 				/* where next would come from */
   for(; s>=cstart; s--, t--) *(t+3)=*(s);
   /* now have fill 0. through 0.00000; use overlaps to avoid tests */
   if (pre>=4) {
-    UINTAT(cstart+pre-4)=UINTAT("0000");
-    UINTAT(cstart)=UINTAT("0.00");
+    memcpy(cstart+pre-4, "0000", 4);
+    memcpy(cstart, "0.00", 4);
     }
    else { /* 2 or 3 */
     *(cstart+pre-1)='0';
-    USHORTAT(cstart)=USHORTAT("0.");
+    memcpy(cstart, "0.", 2);
     }
   *(c+pre)='\0';			/* terminate */
   return string;
@@ -1723,6 +1799,7 @@ decFloat * decFloatZero(decFloat *df){
   void decShowNum(const bcdnum *num, const char *tag) {
     const char *csign="+";		/* sign character */
     uByte *ub;				/* work */
+    uInt  uiwork;			/* for macros */
     if (num->sign==DECFLOAT_Sign) csign="-";
 
     printf(">%s> ", tag);
@@ -1747,7 +1824,7 @@ decFloat * decFloatZero(decFloat *df){
      if (e==0) *c++='0';		/* 0-length case */
       else if (e<1000) {		/* 3 (or fewer) digits case */
        u=&BIN2BCD8[e*4];		/* -> 3 digits + length byte */
-       UINTAT(c)=UINTAT(u+3-*(u+3))|CHARMASK; /* [as above] */
+       UBFROMUI(c, UBTOUI(u+3-*(u+3))|CHARMASK); /* [as above] */
        c+=*(u+3);			/* bump pointer appropriately */
        }
       else {				/* 4-digits */
@@ -1755,7 +1832,7 @@ decFloat * decFloatZero(decFloat *df){
        Int rem=e-(1000*thou);		/* e%1000 */
        *c++=(char)('0'+(char)thou);	/* the thousands digit */
        u=&BIN2BCD8[rem*4];		/* -> 3 digits + length byte */
-       UINTAT(c)=UINTAT(u)|CHARMASK;	/* copy fixed 3+1 characters [is safe] */
+       UBFROMUI(c, UBTOUI(u)|CHARMASK); /* copy fixed 3+1 characters [is safe] */
        c+=3;				/* bump pointer, always 3 digits */
        }
      *c='\0';				/* add terminator */

libdecnumber/decContext.c

@@ -41,12 +41,10 @@
 #include "decContext.h"       /* context and base types */
 #include "decNumberLocal.h"   /* decNumber local types, etc. */
 
-#if DECCHECK
 /* compile-time endian tester [assumes sizeof(Int)>1] */
 static	const  Int mfcone=1;		     /* constant 1 */
-static	const  Flag *mfctop=(Flag *)&mfcone; /* -> top byte */
+static	const  Flag *mfctop=(const Flag *)&mfcone; /* -> top byte */
 #define LITEND *mfctop		   /* named flag; 1=little-endian */
-#endif
 
 /* ------------------------------------------------------------------ */
 /* round-for-reround digits					      */
@@ -80,9 +78,9 @@ decContext *decContextClearStatus(decContext *context, uInt mask) {
 /*  context is the structure to be initialized			      */
 /*  kind selects the required set of default values, one of:	      */
 /*	DEC_INIT_BASE	    -- select ANSI X3-274 defaults	      */
-/*	DEC_INIT_DECIMAL32  -- select IEEE 754r defaults, 32-bit      */
-/*	DEC_INIT_DECIMAL64  -- select IEEE 754r defaults, 64-bit      */
-/*	DEC_INIT_DECIMAL128 -- select IEEE 754r defaults, 128-bit     */
+/*	DEC_INIT_DECIMAL32  -- select IEEE 754 defaults, 32-bit       */
+/*	DEC_INIT_DECIMAL64  -- select IEEE 754 defaults, 64-bit       */
+/*	DEC_INIT_DECIMAL128 -- select IEEE 754 defaults, 128-bit      */
 /*	For any other value a valid context is returned, but with     */
 /*	Invalid_operation set in the status field.		      */
 /*  returns a context structure with the appropriate initial values.  */
@@ -142,15 +140,6 @@ decContext * decContextDefault(decContext *context, Int kind) {
       decContextSetStatus(context, DEC_Invalid_operation); /* trap */
     }
 
-  #if DECCHECK
-  if (LITEND!=DECLITEND) {
-    const char *adj;
-    if (LITEND) adj="little";
-	   else adj="big";
-    printf("Warning: DECLITEND is set to %d, but this computer appears to be %s-endian\n",
-	   DECLITEND, adj);
-    }
-  #endif
   return context;} /* decContextDefault */
 
 /* ------------------------------------------------------------------ */
@@ -390,6 +379,36 @@ const char *decContextStatusToString(const decContext *context) {
   } /* decContextStatusToString */
 
 /* ------------------------------------------------------------------ */
+/* decContextTestEndian -- test whether DECLITEND is set correctly    */
+/*								      */
+/*  quiet is 1 to suppress message; 0 otherwise 		      */
+/*  returns 0 if DECLITEND is correct				      */
+/*	    1 if DECLITEND is incorrect and should be 1 	      */
+/*	   -1 if DECLITEND is incorrect and should be 0 	      */
+/*								      */
+/* A message is displayed if the return value is not 0 and quiet==0.  */
+/*								      */
+/* No error is possible.					      */
+/* ------------------------------------------------------------------ */
+Int decContextTestEndian(Flag quiet) {
+  Int res=0;		      /* optimist */
+  uInt dle=(uInt)DECLITEND;   /* unsign */
+  if (dle>1) dle=1;	      /* ensure 0 or 1 */
+
+  if (LITEND!=DECLITEND) {
+    const char *adj;
+    if (!quiet) {
+      if (LITEND) adj="little";
+	     else adj="big";
+      printf("Warning: DECLITEND is set to %d, but this computer appears to be %s-endian\n",
+	     DECLITEND, adj);
+      }
+    res=(Int)LITEND-dle;
+    }
+  return res;
+  } /* decContextTestEndian */
+
+/* ------------------------------------------------------------------ */
 /* decContextTestSavedStatus -- test bits in saved status	      */
 /*								      */
 /*  oldstatus is the status word to be tested			      */

libdecnumber/decContext.h

@@ -56,15 +56,21 @@
   #define DECCFULLNAME "Decimal Context Descriptor"   /* Verbose name */
   #define DECCAUTHOR   "Mike Cowlishaw" 	      /* Who to blame */
 
-  #include "gstdint.h"		   /* C99 standard integers	      */
+  #if !defined(int32_t)
+    #include <stdint.h> 	   /* C99 standard integers	      */
+  #endif
   #include <stdio.h>		   /* for printf, etc.		      */
   #include <signal.h>		   /* for traps 		      */
 
   /* Extended flags setting -- set this to 0 to use only IEEE flags   */
+  #if !defined(DECEXTFLAG)
   #define DECEXTFLAG 1		   /* 1=enable extended flags	      */
+  #endif
 
   /* Conditional code flag -- set this to 0 for best performance      */
+  #if !defined(DECSUBSET)
   #define DECSUBSET  0		   /* 1=enable subset arithmetic      */
+  #endif
 
   /* Context for operations, with associated constants		      */
   enum rounding {
@@ -102,9 +108,9 @@
   #define DEC_MIN_EMIN	-999999999
   #define DEC_MAX_MATH	    999999 /* max emax, etc., for math funcs. */
 
-  /* Classifications for decimal numbers, aligned with 754r (note     */
-  /* that 'normal' and 'subnormal' are meaningful only with a	      */
-  /* decContext or a fixed size format).			      */
+  /* Classifications for decimal numbers, aligned with 754 (note that */
+  /* 'normal' and 'subnormal' are meaningful only with a decContext   */
+  /* or a fixed size format).					      */
   enum decClass {
     DEC_CLASS_SNAN,
     DEC_CLASS_QNAN,
@@ -170,30 +176,30 @@
     #define DEC_Underflow	     0x00000004
   #endif
 
-  /* IEEE 854 groupings for the flags				      */
+  /* IEEE 754 groupings for the flags				      */
   /* [DEC_Clamped, DEC_Lost_digits, DEC_Rounded, and DEC_Subnormal    */
-  /* are not in IEEE 854]					      */
-  #define DEC_IEEE_854_Division_by_zero	 (DEC_Division_by_zero)
+  /* are not in IEEE 754]					      */
+  #define DEC_IEEE_754_Division_by_zero  (DEC_Division_by_zero)
   #if DECSUBSET
-  #define DEC_IEEE_854_Inexact		 (DEC_Inexact | DEC_Lost_digits)
+  #define DEC_IEEE_754_Inexact		 (DEC_Inexact | DEC_Lost_digits)
   #else
-  #define DEC_IEEE_854_Inexact		 (DEC_Inexact)
+  #define DEC_IEEE_754_Inexact		 (DEC_Inexact)
   #endif
-  #define DEC_IEEE_854_Invalid_operation (DEC_Conversion_syntax |     \
+  #define DEC_IEEE_754_Invalid_operation (DEC_Conversion_syntax |     \
 					  DEC_Division_impossible |   \
 					  DEC_Division_undefined |    \
 					  DEC_Insufficient_storage |  \
 					  DEC_Invalid_context |       \
 					  DEC_Invalid_operation)
-  #define DEC_IEEE_854_Overflow		 (DEC_Overflow)
-  #define DEC_IEEE_854_Underflow	 (DEC_Underflow)
+  #define DEC_IEEE_754_Overflow 	 (DEC_Overflow)
+  #define DEC_IEEE_754_Underflow	 (DEC_Underflow)
 
   /* flags which are normally errors (result is qNaN, infinite, or 0) */
-  #define DEC_Errors (DEC_IEEE_854_Division_by_zero |		      \
-		      DEC_IEEE_854_Invalid_operation |		      \
-		      DEC_IEEE_854_Overflow | DEC_IEEE_854_Underflow)
+  #define DEC_Errors (DEC_IEEE_754_Division_by_zero |		      \
+		      DEC_IEEE_754_Invalid_operation |		      \
+		      DEC_IEEE_754_Overflow | DEC_IEEE_754_Underflow)
   /* flags which cause a result to become qNaN			      */
-  #define DEC_NaNs    DEC_IEEE_854_Invalid_operation
+  #define DEC_NaNs    DEC_IEEE_754_Invalid_operation
 
   /* flags which are normally for information only (finite results)   */
   #if DECSUBSET
@@ -203,6 +209,13 @@
   #define DEC_Information (DEC_Clamped | DEC_Rounded | DEC_Inexact)
   #endif
 
+  /* IEEE 854 names (for compatibility with older decNumber versions) */
+  #define DEC_IEEE_854_Division_by_zero  DEC_IEEE_754_Division_by_zero
+  #define DEC_IEEE_854_Inexact		 DEC_IEEE_754_Inexact
+  #define DEC_IEEE_854_Invalid_operation DEC_IEEE_754_Invalid_operation
+  #define DEC_IEEE_854_Overflow 	 DEC_IEEE_754_Overflow
+  #define DEC_IEEE_854_Underflow	 DEC_IEEE_754_Underflow
+
   /* Name strings for the exceptional conditions		      */
   #define DEC_Condition_CS "Conversion syntax"
   #define DEC_Condition_DZ "Division by zero"
@@ -251,6 +264,7 @@
   extern decContext  * decContextSetStatusFromStringQuiet(decContext *, const char *);
   extern decContext  * decContextSetStatusQuiet(decContext *, uint32_t);
   extern const char  * decContextStatusToString(const decContext *);
+  extern int32_t       decContextTestEndian(uint8_t);
   extern uint32_t      decContextTestSavedStatus(uint32_t, uint32_t);
   extern uint32_t      decContextTestStatus(decContext *, uint32_t);
   extern decContext  * decContextZeroStatus(decContext *);

libdecnumber/decDPD.h

@@ -30,10 +30,9 @@
 
 /* ------------------------------------------------------------------------ */
 /* Binary Coded Decimal and Densely Packed Decimal conversion lookup tables */
-/* [Automatically generated -- do not edit.  2007.05.05]		    */
+/* [Automatically generated -- do not edit.  2008.06.21]		    */
 /* ------------------------------------------------------------------------ */
-/* ------------------------------------------------------------------------ */
-/* For details, see: http://www2.hursley.ibm.com/decimal/DPDecimal.html	    */
+/* For details, see DPDecimal.html on the General Decimal Arithmetic page.  */
 
 #include "decDPDSymbols.h"
 

libdecnumber/decDouble.c

@@ -62,6 +62,7 @@
 #define decFloatFromBCD 	   decDoubleFromBCD
 #define decFloatFromInt32	   decDoubleFromInt32
 #define decFloatFromPacked	   decDoubleFromPacked
+#define decFloatFromPackedChecked  decDoubleFromPackedChecked
 #define decFloatFromString	   decDoubleFromString
 #define decFloatFromUInt32	   decDoubleFromUInt32
 #define decFloatFromWider	   decDoubleFromWider
@@ -145,10 +146,7 @@
 #define decFloatSameQuantum	   decDoubleSameQuantum
 #define decFloatVersion 	   decDoubleVersion
 
-
 #include "decNumberLocal.h"   /* local includes (need DECPMAX) */
 #include "decCommon.c"	      /* non-arithmetic decFloat routines */
 #include "decBasic.c"	      /* basic formats routines */
 
-/* Below here will move to shared file as completed */
-

libdecnumber/decDouble.h

@@ -31,8 +31,6 @@
 /* ------------------------------------------------------------------ */
 /* decDouble.h -- Decimal 64-bit format module header		      */
 /* ------------------------------------------------------------------ */
-/* Please see decFloats.h for an overview and documentation details.  */
-/* ------------------------------------------------------------------ */
 
 #if !defined(DECDOUBLE)
   #define DECDOUBLE
@@ -58,11 +56,14 @@
   #include "decContext.h"
   #include "decQuad.h"
 
-  /* The decDouble decimal 64-bit type, accessible by various types */
+  /* The decDouble decimal 64-bit type, accessible by all sizes */
   typedef union {
     uint8_t   bytes[DECDOUBLE_Bytes];	/* fields: 1, 5, 8, 50 bits */
     uint16_t shorts[DECDOUBLE_Bytes/2];
     uint32_t  words[DECDOUBLE_Bytes/4];
+    #if DECUSE64
+    uint64_t  longs[DECDOUBLE_Bytes/8];
+    #endif
     } decDouble;
 
   /* ---------------------------------------------------------------- */
@@ -75,6 +76,7 @@
   extern decDouble * decDoubleFromBCD(decDouble *, int32_t, const uint8_t *, int32_t);
   extern decDouble * decDoubleFromInt32(decDouble *, int32_t);
   extern decDouble * decDoubleFromPacked(decDouble *, int32_t, const uint8_t *);
+  extern decDouble * decDoubleFromPackedChecked(decDouble *, int32_t, const uint8_t *);
   extern decDouble * decDoubleFromString(decDouble *, const char *, decContext *);
   extern decDouble * decDoubleFromUInt32(decDouble *, uint32_t);
   extern decDouble * decDoubleFromWider(decDouble *, const decQuad *, decContext *);
@@ -160,7 +162,8 @@
 
   /* decNumber conversions; these are implemented as macros so as not  */
   /* to force a dependency on decimal64 and decNumber in decDouble.    */
+  /* decDoubleFromNumber returns a decimal64 * to avoid warnings.      */
   #define decDoubleToNumber(dq, dn) decimal64ToNumber((decimal64 *)(dq), dn)
-  #define decDoubleFromNumber(dq, dn, set) (decDouble *)decimal64FromNumber((decimal64 *)(dq), dn, set)
+  #define decDoubleFromNumber(dq, dn, set) decimal64FromNumber((decimal64 *)(dq), dn, set)
 
 #endif

libdecnumber/decNumber.c

@@ -31,9 +31,9 @@
 /* ------------------------------------------------------------------ */
 /* Decimal Number arithmetic module				      */
 /* ------------------------------------------------------------------ */
-/* This module comprises the routines for General Decimal Arithmetic  */
-/* as defined in the specification which may be found on the	      */
-/* http://www2.hursley.ibm.com/decimal web pages.  It implements both */
+/* This module comprises the routines for arbitrary-precision General */
+/* Decimal Arithmetic as defined in the specification which may be    */
+/* found on the General Decimal Arithmetic pages.  It implements both */
 /* the full ('extended') arithmetic and the simpler ('subset')	      */
 /* arithmetic.							      */
 /*								      */
@@ -41,10 +41,21 @@
 /*								      */
 /* 1. This code is ANSI C89 except:				      */
 /*								      */
-/*       If DECDPUN>4 or DECUSE64=1, the C99 64-bit int64_t and	      */
+/*    a) C99 line comments (double forward slash) are used.  (Most C  */
+/*	 compilers accept these.  If yours does not, a simple script  */
+/*	 can be used to convert them to ANSI C comments.)	      */
+/*								      */
+/*    b) Types from C99 stdint.h are used.  If you do not have this   */
+/*	 header file, see the User's Guide section of the decNumber   */
+/*	 documentation; this lists the necessary definitions.	      */
+/*								      */
+/*    c) If DECDPUN>4 or DECUSE64=1, the C99 64-bit int64_t and       */
 /*	 uint64_t types may be used.  To avoid these, set DECUSE64=0  */
 /*	 and DECDPUN<=4 (see documentation).			      */
 /*								      */
+/*    The code also conforms to C99 restrictions; in particular,      */
+/*    strict aliasing rules are observed.			      */
+/*								      */
 /* 2. The decNumber format which this library uses is optimized for   */
 /*    efficient processing of relatively short numbers; in particular */
 /*    it allows the use of fixed sized structures and minimizes copy  */
@@ -262,7 +273,7 @@ static Int	   decShiftToLeast(Unit *, Int, Int);
 static Int	   decShiftToMost(Unit *, Int, Int);
 static void	   decStatus(decNumber *, uInt, decContext *);
 static void	   decToString(const decNumber *, char[], Flag);
-static decNumber * decTrim(decNumber *, decContext *, Flag, Int *);
+static decNumber * decTrim(decNumber *, decContext *, Flag, Flag, Int *);
 static Int	   decUnitAddSub(const Unit *, Int, const Unit *, Int, Int,
 			      Unit *, Int);
 static Int	   decUnitCompare(const Unit *, Int, const Unit *, Int, Int);
@@ -1156,7 +1167,7 @@ decNumber * decNumberFMA(decNumber *res, const decNumber *lhs,
       }
     #if DECCHECK
      else { /* multiply was OK */
-      if (status!=0) printf("Status=%08lx after FMA multiply\n", status);
+      if (status!=0) printf("Status=%08lx after FMA multiply\n", (LI)status);
       }
     #endif
     /* add the third operand and result -> res, and all is done */
@@ -1301,7 +1312,7 @@ decNumber * decNumberLn(decNumber *res, const decNumber *rhs,
   } /* decNumberLn */
 
 /* ------------------------------------------------------------------ */
-/* decNumberLogB - get adjusted exponent, by 754r rules		      */
+/* decNumberLogB - get adjusted exponent, by 754 rules		      */
 /*								      */
 /*   This computes C = adjustedexponent(A)			      */
 /*								      */
@@ -1338,7 +1349,7 @@ decNumber * decNumberLogB(decNumber *res, const decNumber *rhs,
    else if (decNumberIsZero(rhs)) {
     decNumberZero(res); 		/* prepare for Infinity */
     res->bits=DECNEG|DECINF;		/* -Infinity */
-    status|=DEC_Division_by_zero;	/* as per 754r */
+    status|=DEC_Division_by_zero;	/* as per 754 */
     }
    else { /* finite non-zero */
     Int ae=rhs->exponent+rhs->digits-1; /* adjusted exponent */
@@ -1519,7 +1530,7 @@ decNumber * decNumberLog10(decNumber *res, const decNumber *rhs,
 /* ------------------------------------------------------------------ */
 /* decNumberMax -- compare two Numbers and return the maximum	      */
 /*								      */
-/*   This computes C = A ? B, returning the maximum by 754R rules     */
+/*   This computes C = A ? B, returning the maximum by 754 rules      */
 /*								      */
 /*   res is C, the result.  C may be A and/or B (e.g., X=X?X)	      */
 /*   lhs is A							      */
@@ -1542,7 +1553,7 @@ decNumber * decNumberMax(decNumber *res, const decNumber *lhs,
 /* ------------------------------------------------------------------ */
 /* decNumberMaxMag -- compare and return the maximum by magnitude     */
 /*								      */
-/*   This computes C = A ? B, returning the maximum by 754R rules     */
+/*   This computes C = A ? B, returning the maximum by 754 rules      */
 /*								      */
 /*   res is C, the result.  C may be A and/or B (e.g., X=X?X)	      */
 /*   lhs is A							      */
@@ -1565,7 +1576,7 @@ decNumber * decNumberMaxMag(decNumber *res, const decNumber *lhs,
 /* ------------------------------------------------------------------ */
 /* decNumberMin -- compare two Numbers and return the minimum	      */
 /*								      */
-/*   This computes C = A ? B, returning the minimum by 754R rules     */
+/*   This computes C = A ? B, returning the minimum by 754 rules      */
 /*								      */
 /*   res is C, the result.  C may be A and/or B (e.g., X=X?X)	      */
 /*   lhs is A							      */
@@ -1588,7 +1599,7 @@ decNumber * decNumberMin(decNumber *res, const decNumber *lhs,
 /* ------------------------------------------------------------------ */
 /* decNumberMinMag -- compare and return the minimum by magnitude     */
 /*								      */
-/*   This computes C = A ? B, returning the minimum by 754R rules     */
+/*   This computes C = A ? B, returning the minimum by 754 rules      */
 /*								      */
 /*   res is C, the result.  C may be A and/or B (e.g., X=X?X)	      */
 /*   lhs is A							      */
@@ -1650,7 +1661,7 @@ decNumber * decNumberMinus(decNumber *res, const decNumber *rhs,
 /*   rhs is A							      */
 /*   set is the context 					      */
 /*								      */
-/* This is a generalization of 754r NextDown.			      */
+/* This is a generalization of 754 NextDown.			      */
 /* ------------------------------------------------------------------ */
 decNumber * decNumberNextMinus(decNumber *res, const decNumber *rhs,
 			       decContext *set) {
@@ -1686,7 +1697,7 @@ decNumber * decNumberNextMinus(decNumber *res, const decNumber *rhs,
 /*   rhs is A							      */
 /*   set is the context 					      */
 /*								      */
-/* This is a generalization of 754r NextUp.			      */
+/* This is a generalization of 754 NextUp.			      */
 /* ------------------------------------------------------------------ */
 decNumber * decNumberNextPlus(decNumber *res, const decNumber *rhs,
 			      decContext *set) {
@@ -1718,14 +1729,15 @@ decNumber * decNumberNextPlus(decNumber *res, const decNumber *rhs,
 /* decNumberNextToward -- next towards rhs			      */
 /*								      */
 /*   This computes C = A +/- infinitesimal, rounded towards	      */
-/*   +/-Infinity in the direction of B, as per 754r nextafter rules   */
+/*   +/-Infinity in the direction of B, as per 754-1985 nextafter     */
+/*   modified during revision but dropped from 754-2008.	      */
 /*								      */
 /*   res is C, the result.  C may be A or B.			      */
 /*   lhs is A							      */
 /*   rhs is B							      */
 /*   set is the context 					      */
 /*								      */
-/* This is a generalization of 754r NextAfter.			      */
+/* This is a generalization of 754-1985 NextAfter.		      */
 /* ------------------------------------------------------------------ */
 decNumber * decNumberNextToward(decNumber *res, const decNumber *lhs,
 				const decNumber *rhs, decContext *set) {
@@ -2242,7 +2254,7 @@ decNumber * decNumberPower(decNumber *res, const decNumber *lhs,
     decCopyFit(res, dac, set, &residue, &status);
     decFinish(res, set, &residue, &status);  /* final cleanup */
     #if DECSUBSET
-    if (!set->extended) decTrim(res, set, 0, &dropped); /* trailing zeros */
+    if (!set->extended) decTrim(res, set, 0, 1, &dropped); /* trailing zeros */
     #endif
     } while(0); 			/* end protected */
 
@@ -2337,7 +2349,8 @@ decNumber * decNumberReduce(decNumber *res, const decNumber *rhs,
     /* reduce result to the requested length and copy to result */
     decCopyFit(res, rhs, set, &residue, &status); /* copy & round */
     decFinish(res, set, &residue, &status);	  /* cleanup/set flags */
-    decTrim(res, set, 1, &dropped);		  /* normalize in place */
+    decTrim(res, set, 1, 0, &dropped);		  /* normalize in place */
+						  /* [may clamp] */
     } while(0); 			     /* end protected */
 
   #if DECSUBSET
@@ -2849,7 +2862,7 @@ decNumber * decNumberSquareRoot(decNumber *res, const decNumber *rhs,
       }
 
     /* calculate the ideal (preferred) exponent [floor(exp/2)] */
-    /* [We would like to write: ideal=rhs->exponent>>1, but this */
+    /* [It would be nicer to write: ideal=rhs->exponent>>1, but this */
     /* generates a compiler warning.  Generated code is the same.] */
     ideal=(rhs->exponent&~1)/2; 	/* target */
 
@@ -2876,6 +2889,7 @@ decNumber * decNumberSquareRoot(decNumber *res, const decNumber *rhs,
     /*	 b -- intermediate temporary result (same size as a) */
     /* if any is too long for local storage, then allocate */
     workp=MAXI(set->digits+1, rhs->digits);  /* actual rounding precision */
+    workp=MAXI(workp, 7);		     /* at least 7 for low cases */
     maxp=workp+2;			     /* largest working precision */
 
     needbytes=sizeof(decNumber)+(D2U(rhs->digits)-1)*sizeof(Unit);
@@ -2905,6 +2919,8 @@ decNumber * decNumberSquareRoot(decNumber *res, const decNumber *rhs,
 
     /* set up working context */
     decContextDefault(&workset, DEC_INIT_DECIMAL64);
+    workset.emax=DEC_MAX_EMAX;
+    workset.emin=DEC_MIN_EMIN;
 
     /* [Until further notice, no error is possible and status bits */
     /* (Rounded, etc.) should be ignored, not accumulated.] */
@@ -2945,6 +2961,7 @@ decNumber * decNumberSquareRoot(decNumber *res, const decNumber *rhs,
 	a->lsu[0]=9; a->lsu[1]=5; a->lsu[2]=2;
       #endif
       }
+
     decMultiplyOp(a, a, f, &workset, &ignore);	  /* a=a*f */
     decAddOp(a, a, t, &workset, 0, &ignore);	  /* ..+t */
     /* [a is now the initial approximation for sqrt(f), calculated with */
@@ -2956,16 +2973,14 @@ decNumber * decNumberSquareRoot(decNumber *res, const decNumber *rhs,
     t->lsu[0]=5;			     /* .. */
     t->exponent=-1;			     /* .. */
     workset.digits=3;			     /* initial p */
-    for (;;) {
+    for (; workset.digits<maxp;) {
       /* set p to min(2*p - 2, maxp)  [hence 3; or: 4, 6, 10, ... , maxp] */
-      workset.digits=workset.digits*2-2;
-      if (workset.digits>maxp) workset.digits=maxp;
+      workset.digits=MINI(workset.digits*2-2, maxp);
       /* a = 0.5 * (a + f/a) */
       /* [calculated at p then rounded to currentprecision] */
       decDivideOp(b, f, a, &workset, DIVIDE, &ignore); /* b=f/a */
       decAddOp(b, b, a, &workset, 0, &ignore);	       /* b=b+a */
       decMultiplyOp(a, b, t, &workset, &ignore);       /* a=b*0.5 */
-      if (a->digits==maxp) break;	     /* have required digits */
       } /* loop */
 
     /* Here, 0.1 <= a < 1 [Hull], and a has maxp digits */
@@ -2975,7 +2990,6 @@ decNumber * decNumberSquareRoot(decNumber *res, const decNumber *rhs,
     approxset=*set;			     /* get emin, emax, etc. */
     approxset.round=DEC_ROUND_HALF_EVEN;
     a->exponent+=exp/2; 		     /* set correct exponent */
-
     rstatus=0;				     /* clear status */
     residue=0;				     /* .. and accumulator */
     decCopyFit(a, a, &approxset, &residue, &rstatus);  /* reduce (if needed) */
@@ -3040,13 +3054,13 @@ decNumber * decNumberSquareRoot(decNumber *res, const decNumber *rhs,
     /* count droppable zeros [after any subnormal rounding] by */
     /* trimming a copy */
     decNumberCopy(b, a);
-    decTrim(b, set, 1, &dropped);	     /* [drops trailing zeros] */
+    decTrim(b, set, 1, 1, &dropped);	     /* [drops trailing zeros] */
 
     /* Set Inexact and Rounded.  The answer can only be exact if */
-    /* it is short enough so that squaring it could fit in workp digits, */
-    /* and it cannot have trailing zeros due to clamping, so these are */
-    /* the only (relatively rare) conditions a careful check is needed */
-    if (b->digits*2-1 > workp && !set->clamp) { /* cannot fit */
+    /* it is short enough so that squaring it could fit in workp */
+    /* digits, so this is the only (relatively rare) condition that */
+    /* a careful check is needed */
+    if (b->digits*2-1 > workp) {	     /* cannot fit */
       status|=DEC_Inexact|DEC_Rounded;
       }
      else {				     /* could be exact/unrounded */
@@ -3064,6 +3078,13 @@ decNumber * decNumberSquareRoot(decNumber *res, const decNumber *rhs,
 	  Int todrop=ideal-a->exponent;      /* most that can be dropped */
 	  if (todrop<0) status|=DEC_Rounded; /* ideally would add 0s */
 	   else {			     /* unrounded */
+	    /* there are some to drop, but emax may not allow all */
+	    Int maxexp=set->emax-set->digits+1;
+	    Int maxdrop=maxexp-a->exponent;
+	    if (todrop>maxdrop && set->clamp) { /* apply clamping */
+	      todrop=maxdrop;
+	      status|=DEC_Clamped;
+	      }
 	    if (dropped<todrop) {	     /* clamp to those available */
 	      todrop=dropped;
 	      status|=DEC_Clamped;
@@ -3432,7 +3453,7 @@ decNumber * decNumberCopySign(decNumber *res, const decNumber *lhs,
 /* bcd must have at least dn->digits bytes.  No error is possible; if */
 /* dn is a NaN or Infinite, digits must be 1 and the coefficient 0.   */
 /* ------------------------------------------------------------------ */
-uByte * decNumberGetBCD(const decNumber *dn, uint8_t *bcd) {
+uByte * decNumberGetBCD(const decNumber *dn, uByte *bcd) {
   uByte *ub=bcd+dn->digits-1;	   /* -> lsd */
   const Unit *up=dn->lsu;	   /* Unit pointer, -> lsu */
 
@@ -3532,7 +3553,8 @@ Int decNumberIsSubnormal(const decNumber *dn, decContext *set) {
 /*   returns dn 						      */
 /*								      */
 /* All fields are updated as required.	This is a utility operation,  */
-/* so special values are unchanged and no error is possible.	      */
+/* so special values are unchanged and no error is possible.  The     */
+/* zeros are removed unconditionally.				      */
 /* ------------------------------------------------------------------ */
 decNumber * decNumberTrim(decNumber *dn) {
   Int  dropped; 		   /* work */
@@ -3541,7 +3563,7 @@ decNumber * decNumberTrim(decNumber *dn) {
   if (decCheckOperands(DECUNRESU, DECUNUSED, dn, DECUNCONT)) return dn;
   #endif
   decContextDefault(&set, DEC_INIT_BASE);    /* clamp=0 */
-  return decTrim(dn, &set, 0, &dropped);
+  return decTrim(dn, &set, 0, 1, &dropped);
   } /* decNumberTrim */
 
 /* ------------------------------------------------------------------ */
@@ -3857,7 +3879,7 @@ static decNumber * decAddOp(decNumber *res, const decNumber *lhs,
       #endif
 	/* exponent will be the lower of the two */
 	adjust=lexp-res->exponent;	/* adjustment needed [if -ve] */
-	if (ISZERO(res)) {		/* both 0: special IEEE 854 rules */
+	if (ISZERO(res)) {		/* both 0: special IEEE 754 rules */
 	  if (adjust<0) res->exponent=lexp;  /* set exponent */
 	  /* 0-0 gives +0 unless rounding to -infinity, and -0-0 gives -0 */
 	  if (diffsign) {
@@ -4206,7 +4228,7 @@ static decNumber * decDivideOp(decNumber *res,
   Int	accunits;		   /* count of units accumulated */
   Int	accdigits;		   /* count of digits accumulated */
 
-  Unit	varbuff[SD2U(DECBUFFER*2+DECDPUN)*sizeof(Unit)]; /* buffer for var1 */
+  Unit	varbuff[SD2U(DECBUFFER*2+DECDPUN)];  /* buffer for var1 */
   Unit	*var1=varbuff;		   /* -> var1 array for long subtraction */
   Unit	*varalloc=NULL; 	   /* -> allocated buffer, iff used */
   Unit	*msu1;			   /* -> msu of var1 */
@@ -4761,7 +4783,7 @@ static decNumber * decDivideOp(decNumber *res,
 
     #if DECSUBSET
     /* If a divide then strip trailing zeros if subset [after round] */
-    if (!set->extended && (op==DIVIDE)) decTrim(res, set, 0, &dropped);
+    if (!set->extended && (op==DIVIDE)) decTrim(res, set, 0, 1, &dropped);
     #endif
     } while(0); 			     /* end protected */
 
@@ -5426,7 +5448,7 @@ decNumber * decExpOp(decNumber *res, const decNumber *rhs,
       /* just a sanity check; comment out test to show always */
       if (iterations>p+3)
 	printf("Exp iterations=%ld, status=%08lx, p=%ld, d=%ld\n",
-	       iterations, *status, p, x->digits);
+	       (LI)iterations, (LI)*status, (LI)p, (LI)x->digits);
       #endif
       } /* h<=8 */
 
@@ -5766,7 +5788,7 @@ decNumber * decLnOp(decNumber *res, const decNumber *rhs,
     /* just a sanity check; remove the test to show always */
     if (iterations>24)
       printf("Ln iterations=%ld, status=%08lx, p=%ld, d=%ld\n",
-	    iterations, *status, p, rhs->digits);
+	    (LI)iterations, (LI)*status, (LI)p, (LI)rhs->digits);
     #endif
 
     /* Copy and round the result to res */
@@ -5937,7 +5959,7 @@ static decNumber * decQuantizeOp(decNumber *res, const decNumber *lhs,
       }
      else {
       decFinalize(res, set, &residue, status);	  /* set subnormal flags */
-      *status&=~DEC_Underflow;		/* suppress Underflow [754r] */
+      *status&=~DEC_Underflow;		/* suppress Underflow [as per 754] */
       }
     } while(0); 			/* end protected */
 
@@ -5958,12 +5980,12 @@ static decNumber * decQuantizeOp(decNumber *res, const decNumber *lhs,
 /*     COMPSIG	  -- as COMPARE except that a quiet NaN raises	      */
 /*		     Invalid operation. 			      */
 /*     COMPMAX	  -- returns the larger of the operands, using the    */
-/*		     754r maxnum operation			      */
+/*		     754 maxnum operation			      */
 /*     COMPMAXMAG -- ditto, comparing absolute values		      */
-/*     COMPMIN	  -- the 754r minnum operation			      */
+/*     COMPMIN	  -- the 754 minnum operation			      */
 /*     COMPMINMAG -- ditto, comparing absolute values		      */
 /*     COMTOTAL   -- returns the signum (as a number) giving the      */
-/*		     result of a comparison using 754r total ordering */
+/*		     result of a comparison using 754 total ordering  */
 /*								      */
 /*   res is C, the result.  C may be A and/or B (e.g., X=X?X)	      */
 /*   lhs is A							      */
@@ -6049,7 +6071,7 @@ decNumber * decCompareOp(decNumber *res, const decNumber *lhs,
 
        else if (merged & DECSNAN);	     /* sNaN -> qNaN */
        else { /* here if MIN or MAX and one or two quiet NaNs */
-	/* min or max -- 754r rules ignore single NaN */
+	/* min or max -- 754 rules ignore single NaN */
 	if (!decNumberIsNaN(lhs) || !decNumberIsNaN(rhs)) {
 	  /* just one NaN; force choice to be the non-NaN operand */
 	  op=COMPMAX;
@@ -6091,7 +6113,7 @@ decNumber * decCompareOp(decNumber *res, const decNumber *lhs,
       /* choose the operand for the result */
       const decNumber *choice;
       if (result==0) { /* operands are numerically equal */
-	/* choose according to sign then exponent (see 754r) */
+	/* choose according to sign then exponent (see 754) */
 	uByte slhs=(lhs->bits & DECNEG);
 	uByte srhs=(rhs->bits & DECNEG);
 	#if DECSUBSET
@@ -6556,6 +6578,7 @@ static Int decUnitAddSub(const Unit *a, Int alength,
 /*   dn is the number to trim or normalize			      */
 /*   set is the context to use to check for clamp		      */
 /*   all is 1 to remove all trailing zeros, 0 for just fraction ones  */
+/*   noclamp is 1 to unconditional (unclamped) trim		      */
 /*   dropped returns the number of discarded trailing zeros	      */
 /*   returns dn 						      */
 /*								      */
@@ -6565,7 +6588,7 @@ static Int decUnitAddSub(const Unit *a, Int alength,
 /* so special values are unchanged and no error is possible.	      */
 /* ------------------------------------------------------------------ */
 static decNumber * decTrim(decNumber *dn, decContext *set, Flag all,
-			   Int *dropped) {
+			   Flag noclamp, Int *dropped) {
   Int	d, exp; 		   /* work */
   uInt	cut;			   /* .. */
   Unit	*up;			   /* -> current Unit */
@@ -6611,7 +6634,7 @@ static decNumber * decTrim(decNumber *dn, decContext *set, Flag all,
   if (d==0) return dn;		   /* none to drop */
 
   /* may need to limit drop if clamping */
-  if (set->clamp) {
+  if (set->clamp && !noclamp) {
     Int maxd=set->emax-set->digits+1-dn->exponent;
     if (maxd<=0) return dn;	   /* nothing possible */
     if (d>maxd) d=maxd;
@@ -7319,7 +7342,7 @@ static void decFinalize(decNumber *dn, decContext *set, Int *residue,
 /*								      */
 /* This sets the sign of a number and sets its value to either	      */
 /* Infinity or the maximum finite value, depending on the sign of     */
-/* dn and the rounding mode, following IEEE 854 rules.		      */
+/* dn and the rounding mode, following IEEE 754 rules.		      */
 /* ------------------------------------------------------------------ */
 static void decSetOverflow(decNumber *dn, decContext *set, uInt *status) {
   Flag needmax=0;		   /* result is maximum finite value */
@@ -7403,7 +7426,6 @@ static void decSetMaxValue(decNumber *dn, decContext *set) {
 /* ------------------------------------------------------------------ */
 static void decSetSubnormal(decNumber *dn, decContext *set, Int *residue,
 			    uInt *status) {
-  Int	     dnexp;	      /* saves original exponent */
   decContext workset;	      /* work */
   Int	     etiny, adjust;   /* .. */
 
@@ -7448,7 +7470,6 @@ static void decSetSubnormal(decNumber *dn, decContext *set, Int *residue,
 
   /* adjust>0, so need to rescale the result so exponent becomes Etiny */
   /* [this code is similar to that in rescale] */
-  dnexp=dn->exponent;			/* save exponent */
   workset=*set; 			/* clone rounding, etc. */
   workset.digits=dn->digits-adjust;	/* set requested length */
   workset.emin-=adjust; 		/* and adjust emin to match */
@@ -7456,7 +7477,7 @@ static void decSetSubnormal(decNumber *dn, decContext *set, Int *residue,
   decSetCoeff(dn, &workset, dn->lsu, dn->digits, residue, status);
   decApplyRound(dn, &workset, *residue, status);
 
-  /* Use 754R/854 default rule: Underflow is set iff Inexact */
+  /* Use 754 default rule: Underflow is set iff Inexact */
   /* [independent of whether trapped] */
   if (*status&DEC_Inexact) *status|=DEC_Underflow;
 
@@ -8071,16 +8092,15 @@ static void decCheckInexact(const decNumber *dn, decContext *set) {
 static void *decMalloc(size_t n) {
   uInt	size=n+12;		   /* true size */
   void	*alloc; 		   /* -> allocated storage */
-  uInt	*j;			   /* work */
-  uByte *b, *b0;		   /* .. */
+  uByte *b, *b0;		   /* work */
+  uInt	uiwork; 		   /* for macros */
 
   alloc=malloc(size);		   /* -> allocated storage */
   if (alloc==NULL) return NULL;    /* out of strorage */
   b0=(uByte *)alloc;		   /* as bytes */
   decAllocBytes+=n;		   /* account for storage */
-  j=(uInt *)alloc;		   /* -> first four bytes */
-  *j=n;				   /* save n */
-  /* printf(" alloc ++ dAB: %ld (%d)\n", decAllocBytes, n); */
+  UBFROMUI(alloc, n);		   /* save n */
+  /* printf(" alloc ++ dAB: %ld (%ld)\n", (LI)decAllocBytes, (LI)n); */
   for (b=b0+4; b<b0+8; b++) *b=DECFENCE;
   for (b=b0+n+8; b<b0+n+12; b++) *b=DECFENCE;
   return b0+8;			   /* -> play area */
@@ -8099,20 +8119,20 @@ static void *decMalloc(size_t n) {
 /* is, offset by 8).						      */
 /* ------------------------------------------------------------------ */
 static void decFree(void *alloc) {
-  uInt	*j, n;			   /* pointer, original length */
+  uInt	n;			   /* original length */
   uByte *b, *b0;		   /* work */
+  uInt	uiwork; 		   /* for macros */
 
   if (alloc==NULL) return;	   /* allowed; it's a nop */
   b0=(uByte *)alloc;		   /* as bytes */
   b0-=8;			   /* -> true start of storage */
-  j=(uInt *)b0;			   /* -> first four bytes */
-  n=*j;				   /* lift */
+  n=UBTOUI(b0); 		   /* lift length */
   for (b=b0+4; b<b0+8; b++) if (*b!=DECFENCE)
     printf("=== Corrupt byte [%02x] at offset %d from %ld ===\n", *b,
-	   b-b0-8, (Int)b0);
+	   b-b0-8, (LI)b0);
   for (b=b0+n+8; b<b0+n+12; b++) if (*b!=DECFENCE)
     printf("=== Corrupt byte [%02x] at offset +%d from %ld, n=%ld ===\n", *b,
-	   b-b0-8, (Int)b0, n);
+	   b-b0-8, (LI)b0, (LI)n);
   free(b0);			   /* drop the storage */
   decAllocBytes-=n;		   /* account for storage */
   /* printf(" free -- dAB: %d (%d)\n", decAllocBytes, -n); */

libdecnumber/decNumberLocal.h

@@ -39,34 +39,40 @@
 
 #if !defined(DECNUMBERLOC)
   #define DECNUMBERLOC
-  #define DECVERSION	"decNumber 3.53" /* Package Version [16 max.] */
+  #define DECVERSION	"decNumber 3.61" /* Package Version [16 max.] */
   #define DECNLAUTHOR	"Mike Cowlishaw"	      /* Who to blame */
 
   #include <stdlib.h>	      /* for abs			      */
   #include <string.h>	      /* for memset, strcpy		      */
-  #include "dconfig.h"	      /* for WORDS_BIGENDIAN                  */
 
   /* Conditional code flag -- set this to match hardware platform     */
-  /* 1=little-endian, 0=big-endian	                              */
-  #if WORDS_BIGENDIAN
-  #define DECLITEND 0
-  #else
-  #define DECLITEND 1
+  #if !defined(DECLITEND)
+  #define DECLITEND 1	      /* 1=little-endian, 0=big-endian	      */
   #endif
 
   /* Conditional code flag -- set this to 1 for best performance      */
+  #if !defined(DECUSE64)
   #define DECUSE64  1	      /* 1=use int64s, 0=int32 & smaller only */
+  #endif
 
   /* Conditional check flags -- set these to 0 for best performance   */
+  #if !defined(DECCHECK)
   #define DECCHECK  0	      /* 1 to enable robust checking	      */
+  #endif
+  #if !defined(DECALLOC)
   #define DECALLOC  0	      /* 1 to enable memory accounting	      */
+  #endif
+  #if !defined(DECTRACE)
   #define DECTRACE  0	      /* 1 to trace certain internals, etc.   */
+  #endif
 
   /* Tuning parameter for decNumber (arbitrary precision) module      */
+  #if !defined(DECBUFFER)
   #define DECBUFFER 36	      /* Size basis for local buffers.	This  */
 			      /* should be a common maximum precision */
 			      /* rounded up to a multiple of 4; must  */
 			      /* be zero or positive.		      */
+  #endif
 
   /* ---------------------------------------------------------------- */
   /* Definitions for all modules (general-purpose)		      */
@@ -90,6 +96,7 @@
   /* Development-use definitions				      */
   typedef long int LI;	      /* for printf arguments only	      */
   #define DECNOINT  0	      /* 1 to check no internal use of 'int'  */
+			      /*   or stdint types		      */
   #if DECNOINT
     /* if these interfere with your C includes, do not set DECNOINT   */
     #define int     ?	      /* enable to ensure that plain C 'int'  */
@@ -100,7 +107,6 @@
   extern const uByte  DECSTICKYTAB[10]; /* re-round digits if sticky  */
   extern const uInt   DECPOWERS[10];	/* powers of ten table	      */
   /* The following are included from decDPD.h			      */
-#include "decDPDSymbols.h"
   extern const uShort DPD2BIN[1024];	/* DPD -> 0-999 	      */
   extern const uShort BIN2DPD[1000];	/* 0-999 -> DPD 	      */
   extern const uInt   DPD2BINK[1024];	/* DPD -> 0-999000	      */
@@ -126,17 +132,27 @@
 
   /* ROUNDUP -- round an integer up to a multiple of n		      */
   #define ROUNDUP(i, n) ((((i)+(n)-1)/n)*n)
+  #define ROUNDUP4(i)	(((i)+3)&~3)	/* special for n=4	      */
 
   /* ROUNDDOWN -- round an integer down to a multiple of n	      */
   #define ROUNDDOWN(i, n) (((i)/n)*n)
   #define ROUNDDOWN4(i)   ((i)&~3)	/* special for n=4	      */
 
-  /* References to multi-byte sequences under different sizes	      */
-  /* Refer to a uInt from four bytes starting at a char* or uByte*,   */
-  /* etc.							      */
-  #define UINTAT(b)   (*((uInt	 *)(b)))
-  #define USHORTAT(b) (*((uShort *)(b)))
-  #define UBYTEAT(b)  (*((uByte	 *)(b)))
+  /* References to multi-byte sequences under different sizes; these  */
+  /* require locally declared variables, but do not violate strict    */
+  /* aliasing or alignment (as did the UINTAT simple cast to uInt).   */
+  /* Variables needed are uswork, uiwork, etc. [so do not use at same */
+  /* level in an expression, e.g., UBTOUI(x)==UBTOUI(y) may fail].    */
+
+  /* Return a uInt, etc., from bytes starting at a char* or uByte*    */
+  #define UBTOUS(b)  (memcpy((void *)&uswork, b, 2), uswork)
+  #define UBTOUI(b)  (memcpy((void *)&uiwork, b, 4), uiwork)
+
+  /* Store a uInt, etc., into bytes starting at a char* or uByte*.    */
+  /* Returns i, evaluated, for convenience; has to use uiwork because */
+  /* i may be an expression.					      */
+  #define UBFROMUS(b, i)  (uswork=(i), memcpy(b, (void *)&uswork, 2), uswork)
+  #define UBFROMUI(b, i)  (uiwork=(i), memcpy(b, (void *)&uiwork, 4), uiwork)
 
   /* X10 and X100 -- multiply integer i by 10 or 100		      */
   /* [shifts are usually faster than multiply; could be conditional]  */
@@ -308,12 +324,12 @@
   #define DECWORDS  (DECBYTES/4)
   #define DECWWORDS (DECWBYTES/4)
   #if DECLITEND
-    #define DFWORD(df, off) ((df)->words[DECWORDS-1-(off)])
     #define DFBYTE(df, off)   ((df)->bytes[DECBYTES-1-(off)])
+    #define DFWORD(df, off)   ((df)->words[DECWORDS-1-(off)])
     #define DFWWORD(dfw, off) ((dfw)->words[DECWWORDS-1-(off)])
   #else
-    #define DFWORD(df, off) ((df)->words[off])
     #define DFBYTE(df, off)   ((df)->bytes[off])
+    #define DFWORD(df, off)   ((df)->words[off])
     #define DFWWORD(dfw, off) ((dfw)->words[off])
   #endif
 
@@ -326,7 +342,6 @@
   #define DFISSNAN(df)	  ((DFWORD(df, 0)&0x7e000000)==0x7e000000)
 
   /* Shared lookup tables					      */
-#include "decCommonSymbols.h"
   extern const uInt   DECCOMBMSD[64];	/* Combination field -> MSD   */
   extern const uInt   DECCOMBFROM[48];	/* exp+msd -> Combination     */
 
@@ -412,21 +427,31 @@
       || ((lo)&(((uInt)0x6e)<<(k)))!=(((uInt)0x6e)<<(k)))
 
     /* Macro to test whether a full-length (length DECPMAX) BCD8      */
-    /* coefficient is zero					      */
-    /* test just the LSWord first, then the remainder		      */
+    /* coefficient, starting at uByte u, is all zeros		      */
+    /* Test just the LSWord first, then the remainder as a sequence   */
+    /* of tests in order to avoid same-level use of UBTOUI	      */
     #if DECPMAX==7
-      #define ISCOEFFZERO(u) (UINTAT((u)+DECPMAX-4)==0		    \
-	&& UINTAT((u)+DECPMAX-7)==0)
+      #define ISCOEFFZERO(u) (					    \
+	   UBTOUI((u)+DECPMAX-4)==0				    \
+	&& UBTOUS((u)+DECPMAX-6)==0				    \
+	&& *(u)==0)
     #elif DECPMAX==16
-      #define ISCOEFFZERO(u) (UINTAT((u)+DECPMAX-4)==0		    \
-	&& (UINTAT((u)+DECPMAX-8)+UINTAT((u)+DECPMAX-12)	    \
-	   +UINTAT((u)+DECPMAX-16))==0)
+      #define ISCOEFFZERO(u) (					    \
+	   UBTOUI((u)+DECPMAX-4)==0				    \
+	&& UBTOUI((u)+DECPMAX-8)==0				    \
+	&& UBTOUI((u)+DECPMAX-12)==0				    \
+	&& UBTOUI(u)==0)
     #elif DECPMAX==34
-      #define ISCOEFFZERO(u) (UINTAT((u)+DECPMAX-4)==0		    \
-	&& (UINTAT((u)+DECPMAX-8) +UINTAT((u)+DECPMAX-12)	    \
-	   +UINTAT((u)+DECPMAX-16)+UINTAT((u)+DECPMAX-20)	    \
-	   +UINTAT((u)+DECPMAX-24)+UINTAT((u)+DECPMAX-28)	    \
-	   +UINTAT((u)+DECPMAX-32)+USHORTAT((u)+DECPMAX-34))==0)
+      #define ISCOEFFZERO(u) (					    \
+	   UBTOUI((u)+DECPMAX-4)==0				    \
+	&& UBTOUI((u)+DECPMAX-8)==0				    \
+	&& UBTOUI((u)+DECPMAX-12)==0				    \
+	&& UBTOUI((u)+DECPMAX-16)==0				    \
+	&& UBTOUI((u)+DECPMAX-20)==0				    \
+	&& UBTOUI((u)+DECPMAX-24)==0				    \
+	&& UBTOUI((u)+DECPMAX-28)==0				    \
+	&& UBTOUI((u)+DECPMAX-32)==0				    \
+	&& UBTOUS(u)==0)
     #endif
 
     /* Macros and masks for the exponent continuation field and MSD   */
@@ -448,20 +473,15 @@
     #define ECONNANMASK ((0x01ffffff>>(32-6-DECECONL))<<(32-6-DECECONL))
 
     /* Macros to decode the coefficient in a finite decFloat *df into */
-    /* a BCD string (uByte *bcdin) of length DECPMAX uBytes	      */
+    /* a BCD string (uByte *bcdin) of length DECPMAX uBytes.	      */
 
-    /* In-line sequence to convert 10 bits at right end of uInt dpd   */
-    /* to three BCD8 digits starting at uByte u.  Note that an extra  */
-    /* byte is written to the right of the three digits because this  */
-    /* moves four at a time for speed; the alternative macro moves    */
-    /* exactly three bytes					      */
-    #define dpd2bcd8(u, dpd) {				 \
-      UINTAT(u)=UINTAT(&DPD2BCD8[((dpd)&0x3ff)*4]);}
-
-    #define dpd2bcd83(u, dpd) {				 \
-      *(u)=DPD2BCD8[((dpd)&0x3ff)*4];			 \
-      *(u+1)=DPD2BCD8[((dpd)&0x3ff)*4+1];		 \
-      *(u+2)=DPD2BCD8[((dpd)&0x3ff)*4+2];}
+    /* In-line sequence to convert least significant 10 bits of uInt  */
+    /* dpd to three BCD8 digits starting at uByte u.  Note that an    */
+    /* extra byte is written to the right of the three digits because */
+    /* four bytes are moved at a time for speed; the alternative      */
+    /* macro moves exactly three bytes (usually slower).	      */
+    #define dpd2bcd8(u, dpd)  memcpy(u, &DPD2BCD8[((dpd)&0x3ff)*4], 4)
+    #define dpd2bcd83(u, dpd) memcpy(u, &DPD2BCD8[((dpd)&0x3ff)*4], 3)
 
     /* Decode the declets.  After extracting each one, it is decoded  */
     /* to BCD8 using a table lookup (also used for variable-length    */
@@ -597,8 +617,8 @@
     #endif
 
     /* Macros to decode the coefficient in a finite decFloat *df into */
-    /* a base-thousand uInt array, with the least-significant 0-999   */
-    /* 'digit' at offset 0.					      */
+    /* a base-thousand uInt array (of size DECLETS+1, to allow for    */
+    /* the MSD), with the least-significant 0-999 'digit' at offset 0.*/
 
     /* Decode the declets.  After extracting each one, it is decoded  */
     /* to binary using a table lookup.				      */
@@ -640,9 +660,72 @@
       (buf)[9]=DPD2BIN[((sourhi<<6) | (sourmh>>26))&0x3ff];   \
       (buf)[10]=DPD2BIN[(sourhi>>4)&0x3ff];		      \
       (buf)[11]=DECCOMBMSD[sourhi>>26];}
+    #endif
+
 
+    /* Macros to decode the coefficient in a finite decFloat *df and  */
+    /* add to a base-thousand uInt array (as for GETCOEFFTHOU).       */
+    /* After the addition then most significant 'digit' in the array  */
+    /* might have a value larger then 10 (with a maximum of 19).      */
+    #if DECPMAX==7
+    #define ADDCOEFFTHOU(df, buf) {			      \
+      uInt sourhi=DFWORD(df, 0);			      \
+      (buf)[0]+=DPD2BIN[sourhi&0x3ff];			      \
+      if (buf[0]>999) {buf[0]-=1000; buf[1]++;} 	      \
+      (buf)[1]+=DPD2BIN[(sourhi>>10)&0x3ff];		      \
+      if (buf[1]>999) {buf[1]-=1000; buf[2]++;} 	      \
+      (buf)[2]+=DECCOMBMSD[sourhi>>26];}
+
+    #elif DECPMAX==16
+    #define ADDCOEFFTHOU(df, buf) {			      \
+      uInt sourhi, sourlo;				      \
+      sourlo=DFWORD(df, 1);				      \
+      (buf)[0]+=DPD2BIN[sourlo&0x3ff];			      \
+      if (buf[0]>999) {buf[0]-=1000; buf[1]++;} 	      \
+      (buf)[1]+=DPD2BIN[(sourlo>>10)&0x3ff];		      \
+      if (buf[1]>999) {buf[1]-=1000; buf[2]++;} 	      \
+      (buf)[2]+=DPD2BIN[(sourlo>>20)&0x3ff];		      \
+      if (buf[2]>999) {buf[2]-=1000; buf[3]++;} 	      \
+      sourhi=DFWORD(df, 0);				      \
+      (buf)[3]+=DPD2BIN[((sourhi<<2) | (sourlo>>30))&0x3ff];  \
+      if (buf[3]>999) {buf[3]-=1000; buf[4]++;} 	      \
+      (buf)[4]+=DPD2BIN[(sourhi>>8)&0x3ff];		      \
+      if (buf[4]>999) {buf[4]-=1000; buf[5]++;} 	      \
+      (buf)[5]+=DECCOMBMSD[sourhi>>26];}
+
+    #elif DECPMAX==34
+    #define ADDCOEFFTHOU(df, buf) {			      \
+      uInt sourhi, sourmh, sourml, sourlo;		      \
+      sourlo=DFWORD(df, 3);				      \
+      (buf)[0]+=DPD2BIN[sourlo&0x3ff];			      \
+      if (buf[0]>999) {buf[0]-=1000; buf[1]++;} 	      \
+      (buf)[1]+=DPD2BIN[(sourlo>>10)&0x3ff];		      \
+      if (buf[1]>999) {buf[1]-=1000; buf[2]++;} 	      \
+      (buf)[2]+=DPD2BIN[(sourlo>>20)&0x3ff];		      \
+      if (buf[2]>999) {buf[2]-=1000; buf[3]++;} 	      \
+      sourml=DFWORD(df, 2);				      \
+      (buf)[3]+=DPD2BIN[((sourml<<2) | (sourlo>>30))&0x3ff];  \
+      if (buf[3]>999) {buf[3]-=1000; buf[4]++;} 	      \
+      (buf)[4]+=DPD2BIN[(sourml>>8)&0x3ff];		      \
+      if (buf[4]>999) {buf[4]-=1000; buf[5]++;} 	      \
+      (buf)[5]+=DPD2BIN[(sourml>>18)&0x3ff];		      \
+      if (buf[5]>999) {buf[5]-=1000; buf[6]++;} 	      \
+      sourmh=DFWORD(df, 1);				      \
+      (buf)[6]+=DPD2BIN[((sourmh<<4) | (sourml>>28))&0x3ff];  \
+      if (buf[6]>999) {buf[6]-=1000; buf[7]++;} 	      \
+      (buf)[7]+=DPD2BIN[(sourmh>>6)&0x3ff];		      \
+      if (buf[7]>999) {buf[7]-=1000; buf[8]++;} 	      \
+      (buf)[8]+=DPD2BIN[(sourmh>>16)&0x3ff];		      \
+      if (buf[8]>999) {buf[8]-=1000; buf[9]++;} 	      \
+      sourhi=DFWORD(df, 0);				      \
+      (buf)[9]+=DPD2BIN[((sourhi<<6) | (sourmh>>26))&0x3ff];  \
+      if (buf[9]>999) {buf[9]-=1000; buf[10]++;}	      \
+      (buf)[10]+=DPD2BIN[(sourhi>>4)&0x3ff];		      \
+      if (buf[10]>999) {buf[10]-=1000; buf[11]++;}	      \
+      (buf)[11]+=DECCOMBMSD[sourhi>>26];}
     #endif
 
+
     /* Set a decFloat to the maximum positive finite number (Nmax)    */
     #if DECPMAX==7
     #define DFSETNMAX(df)	     \

libdecnumber/decQuad.c

@@ -56,6 +56,7 @@
 #define decFloatFromBCD 	   decQuadFromBCD
 #define decFloatFromInt32	   decQuadFromInt32
 #define decFloatFromPacked	   decQuadFromPacked
+#define decFloatFromPackedChecked  decQuadFromPackedChecked
 #define decFloatFromString	   decQuadFromString
 #define decFloatFromUInt32	   decQuadFromUInt32
 #define decFloatFromWider	   decQuadFromWider
@@ -139,7 +140,6 @@
 #define decFloatSameQuantum	   decQuadSameQuantum
 #define decFloatVersion 	   decQuadVersion
 
-
 #include "decNumberLocal.h"   /* local includes (need DECPMAX) */
 #include "decCommon.c"	      /* non-arithmetic decFloat routines */
 #include "decBasic.c"	      /* basic formats routines */

libdecnumber/decQuad.h

@@ -31,8 +31,6 @@
 /* ------------------------------------------------------------------ */
 /* decQuad.h -- Decimal 128-bit format module header		      */
 /* ------------------------------------------------------------------ */
-/* Please see decFloats.h for an overview and documentation details.  */
-/* ------------------------------------------------------------------ */
 /* This include file is always included by decSingle and decDouble,   */
 /* and therefore also holds useful constants used by all three.       */
 
@@ -59,11 +57,14 @@
   /* Required include						      */
   #include "decContext.h"
 
-  /* The decQuad decimal 128-bit type, accessible by various types */
+  /* The decQuad decimal 128-bit type, accessible by all sizes */
   typedef union {
     uint8_t   bytes[DECQUAD_Bytes];	/* fields: 1, 5, 12, 110 bits */
     uint16_t shorts[DECQUAD_Bytes/2];
     uint32_t  words[DECQUAD_Bytes/4];
+    #if DECUSE64
+    uint64_t  longs[DECQUAD_Bytes/8];
+    #endif
     } decQuad;
 
   /* ---------------------------------------------------------------- */
@@ -99,6 +100,7 @@
   extern decQuad * decQuadFromBCD(decQuad *, int32_t, const uint8_t *, int32_t);
   extern decQuad * decQuadFromInt32(decQuad *, int32_t);
   extern decQuad * decQuadFromPacked(decQuad *, int32_t, const uint8_t *);
+  extern decQuad * decQuadFromPackedChecked(decQuad *, int32_t, const uint8_t *);
   extern decQuad * decQuadFromString(decQuad *, const char *, decContext *);
   extern decQuad * decQuadFromUInt32(decQuad *, uint32_t);
   extern int32_t   decQuadGetCoefficient(const decQuad *, uint8_t *);
@@ -182,7 +184,8 @@
 
   /* decNumber conversions; these are implemented as macros so as not  */
   /* to force a dependency on decimal128 and decNumber in decQuad.     */
+  /* decQuadFromNumber returns a decimal128 * to avoid warnings.       */
   #define decQuadToNumber(dq, dn) decimal128ToNumber((decimal128 *)(dq), dn)
-  #define decQuadFromNumber(dq, dn, set) (decQuad *)decimal128FromNumber((decimal128 *)(dq), dn, set)
+  #define decQuadFromNumber(dq, dn, set) decimal128FromNumber((decimal128 *)(dq), dn, set)
 
 #endif

libdecnumber/decSingle.c

@@ -31,8 +31,6 @@
 /* ------------------------------------------------------------------ */
 /* decSingle.c -- decSingle operations module			      */
 /* ------------------------------------------------------------------ */
-/* This module comprises decSingle operations (including conversions) */
-/* ------------------------------------------------------------------ */
 
 #include "decContext.h"       /* public includes */
 #include "decSingle.h"	      /* public includes */
@@ -61,6 +59,7 @@
 /* Utility (binary results, extractors, etc.) */
 #define decFloatFromBCD 	   decSingleFromBCD
 #define decFloatFromPacked	   decSingleFromPacked
+#define decFloatFromPackedChecked  decSingleFromPackedChecked
 #define decFloatFromString	   decSingleFromString
 #define decFloatFromWider	   decSingleFromWider
 #define decFloatGetCoefficient	   decSingleGetCoefficient

libdecnumber/decSingle.h

@@ -31,8 +31,6 @@
 /* ------------------------------------------------------------------ */
 /* decSingle.h -- Decimal 32-bit format module header		      */
 /* ------------------------------------------------------------------ */
-/* Please see decFloats.h for an overview and documentation details.  */
-/* ------------------------------------------------------------------ */
 
 #if !defined(DECSINGLE)
   #define DECSINGLE
@@ -59,7 +57,7 @@
   #include "decQuad.h"
   #include "decDouble.h"
 
-  /* The decSingle decimal 32-bit type, accessible by various types */
+  /* The decSingle decimal 32-bit type, accessible by all sizes */
   typedef union {
     uint8_t   bytes[DECSINGLE_Bytes];	/* fields: 1, 5, 6, 20 bits */
     uint16_t shorts[DECSINGLE_Bytes/2];
@@ -75,6 +73,7 @@
   /* Utilities (binary argument(s) or result, extractors, etc.) */
   extern decSingle * decSingleFromBCD(decSingle *, int32_t, const uint8_t *, int32_t);
   extern decSingle * decSingleFromPacked(decSingle *, int32_t, const uint8_t *);
+  extern decSingle * decSingleFromPackedChecked(decSingle *, int32_t, const uint8_t *);
   extern decSingle * decSingleFromString(decSingle *, const char *, decContext *);
   extern decSingle * decSingleFromWider(decSingle *, const decDouble *, decContext *);
   extern int32_t     decSingleGetCoefficient(const decSingle *, uint8_t *);
@@ -97,7 +96,8 @@
 
   /* decNumber conversions; these are implemented as macros so as not  */
   /* to force a dependency on decimal32 and decNumber in decSingle.    */
+  /* decSingleFromNumber returns a decimal32 * to avoid warnings.      */
   #define decSingleToNumber(dq, dn) decimal32ToNumber((decimal32 *)(dq), dn)
-  #define decSingleFromNumber(dq, dn, set) (decSingle *)decimal32FromNumber((decimal32 *)(dq), dn, set)
+  #define decSingleFromNumber(dq, dn, set) decimal32FromNumber((decimal32 *)(dq), dn, set)
 
 #endif

libdecnumber/dpd/decimal128.c

@@ -89,8 +89,8 @@ decimal128 * decimal128FromNumber(decimal128 *d128, const decNumber *dn,
   Int ae;			   /* adjusted exponent */
   decNumber  dw;		   /* work */
   decContext dc;		   /* .. */
-  uInt *pu;			   /* .. */
   uInt comb, exp;		   /* .. */
+  uInt uiwork;			   /* for macros */
   uInt targar[4]={0,0,0,0};	   /* target 128-bit */
   #define targhi targar[3]	   /* name the word with the sign */
   #define targmh targar[2]	   /* name the words */
@@ -172,18 +172,19 @@ decimal128 * decimal128FromNumber(decimal128 *d128, const decNumber *dn,
   if (dn->bits&DECNEG) targhi|=0x80000000; /* add sign bit */
 
   /* now write to storage; this is endian */
-  pu=(uInt *)d128->bytes;	   /* overlay */
   if (DECLITEND) {
-    pu[0]=targlo;		   /* directly store the low int */
-    pu[1]=targml;		   /* then the mid-low */
-    pu[2]=targmh;		   /* then the mid-high */
-    pu[3]=targhi;		   /* then the high int */
+    /* lo -> hi */
+    UBFROMUI(d128->bytes,    targlo);
+    UBFROMUI(d128->bytes+4,  targml);
+    UBFROMUI(d128->bytes+8,  targmh);
+    UBFROMUI(d128->bytes+12, targhi);
     }
    else {
-    pu[0]=targhi;		   /* directly store the high int */
-    pu[1]=targmh;		   /* then the mid-high */
-    pu[2]=targml;		   /* then the mid-low */
-    pu[3]=targlo;		   /* then the low int */
+    /* hi -> lo */
+    UBFROMUI(d128->bytes,    targhi);
+    UBFROMUI(d128->bytes+4,  targmh);
+    UBFROMUI(d128->bytes+8,  targml);
+    UBFROMUI(d128->bytes+12, targlo);
     }
 
   if (status!=0) decContextSetStatus(set, status); /* pass on status */
@@ -201,8 +202,8 @@ decNumber * decimal128ToNumber(const decimal128 *d128, decNumber *dn) {
   uInt msd;			   /* coefficient MSD */
   uInt exp;			   /* exponent top two bits */
   uInt comb;			   /* combination field */
-  const uInt *pu;		   /* work */
-  Int  need;			   /* .. */
+  Int  need;			   /* work */
+  uInt uiwork;			   /* for macros */
   uInt sourar[4];		   /* source 128-bit */
   #define sourhi sourar[3]	   /* name the word with the sign */
   #define sourmh sourar[2]	   /* and the mid-high word */
@@ -210,18 +211,17 @@ decNumber * decimal128ToNumber(const decimal128 *d128, decNumber *dn) {
   #define sourlo sourar[0]	   /* and the lowest word */
 
   /* load source from storage; this is endian */
-  pu=(const uInt *)d128->bytes;	   /* overlay */
   if (DECLITEND) {
-    sourlo=pu[0];		   /* directly load the low int */
-    sourml=pu[1];		   /* then the mid-low */
-    sourmh=pu[2];		   /* then the mid-high */
-    sourhi=pu[3];		   /* then the high int */
+    sourlo=UBTOUI(d128->bytes	); /* directly load the low int */
+    sourml=UBTOUI(d128->bytes+4 ); /* then the mid-low */
+    sourmh=UBTOUI(d128->bytes+8 ); /* then the mid-high */
+    sourhi=UBTOUI(d128->bytes+12); /* then the high int */
     }
    else {
-    sourhi=pu[0];		   /* directly load the high int */
-    sourmh=pu[1];		   /* then the mid-high */
-    sourml=pu[2];		   /* then the mid-low */
-    sourlo=pu[3];		   /* then the low int */
+    sourhi=UBTOUI(d128->bytes	); /* directly load the high int */
+    sourmh=UBTOUI(d128->bytes+4 ); /* then the mid-high */
+    sourml=UBTOUI(d128->bytes+8 ); /* then the mid-low */
+    sourlo=UBTOUI(d128->bytes+12); /* then the low int */
     }
 
   comb=(sourhi>>26)&0x1f;	   /* combination field */
@@ -291,11 +291,11 @@ char * decimal128ToString(const decimal128 *d128, char *string){
   uInt comb;			   /* combination field */
   char *cstart; 		   /* coefficient start */
   char *c;			   /* output pointer in string */
-  const uInt *pu;		   /* work */
+  const uByte *u;		   /* work */
   char *s, *t;			   /* .. (source, target) */
   Int  dpd;			   /* .. */
   Int  pre, e;			   /* .. */
-  const uByte *u;		   /* .. */
+  uInt uiwork;			   /* for macros */
 
   uInt sourar[4];		   /* source 128-bit */
   #define sourhi sourar[3]	   /* name the word with the sign */
@@ -304,18 +304,17 @@ char * decimal128ToString(const decimal128 *d128, char *string){
   #define sourlo sourar[0]	   /* and the lowest word */
 
   /* load source from storage; this is endian */
-  pu=(const uInt *)d128->bytes;	   /* overlay */
   if (DECLITEND) {
-    sourlo=pu[0];		   /* directly load the low int */
-    sourml=pu[1];		   /* then the mid-low */
-    sourmh=pu[2];		   /* then the mid-high */
-    sourhi=pu[3];		   /* then the high int */
+    sourlo=UBTOUI(d128->bytes	); /* directly load the low int */
+    sourml=UBTOUI(d128->bytes+4 ); /* then the mid-low */
+    sourmh=UBTOUI(d128->bytes+8 ); /* then the mid-high */
+    sourhi=UBTOUI(d128->bytes+12); /* then the high int */
     }
    else {
-    sourhi=pu[0];		   /* directly load the high int */
-    sourmh=pu[1];		   /* then the mid-high */
-    sourml=pu[2];		   /* then the mid-low */
-    sourlo=pu[3];		   /* then the low int */
+    sourhi=UBTOUI(d128->bytes	); /* directly load the high int */
+    sourmh=UBTOUI(d128->bytes+4 ); /* then the mid-high */
+    sourml=UBTOUI(d128->bytes+8 ); /* then the mid-low */
+    sourlo=UBTOUI(d128->bytes+12); /* then the low int */
     }
 
   c=string;			   /* where result will go */
@@ -483,7 +482,7 @@ decimal128 * decimal128FromString(decimal128 *result, const char *string,
 /*   returns 1 if the encoding of d128 is canonical, 0 otherwise      */
 /* No error is possible.					      */
 /* ------------------------------------------------------------------ */
-uint32_t decimal128IsCanonical(const decimal128 *d128) {
+uInt decimal128IsCanonical(const decimal128 *d128) {
   decNumber dn; 			/* work */
   decimal128 canon;			 /* .. */
   decContext dc;			/* .. */
@@ -532,9 +531,9 @@ decimal128 * decimal128Canonical(decimal128 *result, const decimal128 *d128) {
 /* This assumes range has been checked and exponent previously 0; */
 /* type of exponent must be unsigned */
 #define decimal128SetExpCon(d, e) {				      \
-  (d)->bytes[0]|=(uint8_t)((e)>>10);				      \
-  (d)->bytes[1] =(uint8_t)(((e)&0x3fc)>>2);			      \
-  (d)->bytes[2]|=(uint8_t)(((e)&0x03)<<6);}
+  (d)->bytes[0]|=(uByte)((e)>>10);				      \
+  (d)->bytes[1] =(uByte)(((e)&0x3fc)>>2);			      \
+  (d)->bytes[2]|=(uByte)(((e)&0x03)<<6);}
 
 /* ------------------------------------------------------------------ */
 /* decimal128Show -- display a decimal128 in hexadecimal [debug aid]  */

libdecnumber/dpd/decimal128.h

@@ -76,13 +76,13 @@
     #define DECIMAL_Inf     0x78	/* 0 11110 00 Infinity	      */
   #endif
 
-  #include "decimal128Local.h"
+#include "decimal128Local.h"
 
   /* ---------------------------------------------------------------- */
   /* Routines							      */
   /* ---------------------------------------------------------------- */
 
-  #include "decimal128Symbols.h"
+#include "decimal128Symbols.h"
 
   /* String conversions 					      */
   decimal128 * decimal128FromString(decimal128 *, const char *, decContext *);

libdecnumber/dpd/decimal32.c

@@ -89,8 +89,8 @@ decimal32 * decimal32FromNumber(decimal32 *d32, const decNumber *dn,
   Int ae;			   /* adjusted exponent */
   decNumber  dw;		   /* work */
   decContext dc;		   /* .. */
-  uInt *pu;			   /* .. */
   uInt comb, exp;		   /* .. */
+  uInt uiwork;			   /* for macros */
   uInt targ=0;			   /* target 32-bit */
 
   /* If the number has too many digits, or the exponent could be */
@@ -175,8 +175,7 @@ decimal32 * decimal32FromNumber(decimal32 *d32, const decNumber *dn,
   if (dn->bits&DECNEG) targ|=0x80000000;  /* add sign bit */
 
   /* now write to storage; this is endian */
-  pu=(uInt *)d32->bytes;	   /* overlay */
-  *pu=targ;			   /* directly store the int */
+  UBFROMUI(d32->bytes, targ);	   /* directly store the int */
 
   if (status!=0) decContextSetStatus(set, status); /* pass on status */
   /* decimal32Show(d32); */
@@ -194,11 +193,10 @@ decNumber * decimal32ToNumber(const decimal32 *d32, decNumber *dn) {
   uInt exp;			   /* exponent top two bits */
   uInt comb;			   /* combination field */
   uInt sour;			   /* source 32-bit */
-  const uInt *pu;		   /* work */
+  uInt uiwork;			   /* for macros */
 
   /* load source from storage; this is endian */
-  pu=(const uInt *)d32->bytes;	   /* overlay */
-  sour=*pu;			   /* directly load the int */
+  sour=UBTOUI(d32->bytes);	   /* directly load the int */
 
   comb=(sour>>26)&0x1f; 	   /* combination field */
 
@@ -264,16 +262,15 @@ char * decimal32ToString(const decimal32 *d32, char *string){
   uInt comb;			   /* combination field */
   char *cstart; 		   /* coefficient start */
   char *c;			   /* output pointer in string */
-  const uInt *pu;		   /* work */
-  const uByte *u;		   /* .. */
+  const uByte *u;		   /* work */
   char *s, *t;			   /* .. (source, target) */
   Int  dpd;			   /* .. */
   Int  pre, e;			   /* .. */
+  uInt uiwork;			   /* for macros */
   uInt sour;			   /* source 32-bit */
 
   /* load source from storage; this is endian */
-  pu=(const uInt *)d32->bytes;	   /* overlay */
-  sour=*pu;			   /* directly load the int */
+  sour=UBTOUI(d32->bytes);	   /* directly load the int */
 
   c=string;			   /* where result will go */
   if (((Int)sour)<0) *c++='-';	   /* handle sign */
@@ -412,7 +409,7 @@ decimal32 * decimal32FromString(decimal32 *result, const char *string,
 /*   returns 1 if the encoding of d32 is canonical, 0 otherwise       */
 /* No error is possible.					      */
 /* ------------------------------------------------------------------ */
-uint32_t decimal32IsCanonical(const decimal32 *d32) {
+uInt decimal32IsCanonical(const decimal32 *d32) {
   decNumber dn; 			/* work */
   decimal32 canon;			/* .. */
   decContext dc;			/* .. */
@@ -460,8 +457,8 @@ decimal32 * decimal32Canonical(decimal32 *result, const decimal32 *d32) {
 /* This assumes range has been checked and exponent previously 0; */
 /* type of exponent must be unsigned */
 #define decimal32SetExpCon(d, e) {				      \
-  (d)->bytes[0]|=(uint8_t)((e)>>4);				      \
-  (d)->bytes[1]|=(uint8_t)(((e)&0x0F)<<4);}
+  (d)->bytes[0]|=(uByte)((e)>>4);				      \
+  (d)->bytes[1]|=(uByte)(((e)&0x0F)<<4);}
 
 /* ------------------------------------------------------------------ */
 /* decimal32Show -- display a decimal32 in hexadecimal [debug aid]    */

libdecnumber/dpd/decimal32.h

@@ -80,7 +80,7 @@
   /* Routines							      */
   /* ---------------------------------------------------------------- */
 
-  #include "decimal32Symbols.h"
+#include "decimal32Symbols.h"
 
   /* String conversions 					      */
   decimal32 * decimal32FromString(decimal32 *, const char *, decContext *);

libdecnumber/dpd/decimal64.c

@@ -95,8 +95,8 @@ decimal64 * decimal64FromNumber(decimal64 *d64, const decNumber *dn,
   Int ae;			   /* adjusted exponent */
   decNumber  dw;		   /* work */
   decContext dc;		   /* .. */
-  uInt *pu;			   /* .. */
   uInt comb, exp;		   /* .. */
+  uInt uiwork;			   /* for macros */
   uInt targar[2]={0, 0};	   /* target 64-bit */
   #define targhi targar[1]	   /* name the word with the sign */
   #define targlo targar[0]	   /* and the other */
@@ -193,14 +193,15 @@ decimal64 * decimal64FromNumber(decimal64 *d64, const decNumber *dn,
   if (dn->bits&DECNEG) targhi|=0x80000000; /* add sign bit */
 
   /* now write to storage; this is now always endian */
-  pu=(uInt *)d64->bytes;	   /* overlay */
   if (DECLITEND) {
-    pu[0]=targar[0];		   /* directly store the low int */
-    pu[1]=targar[1];		   /* then the high int */
+    /* lo int then hi */
+    UBFROMUI(d64->bytes,   targar[0]);
+    UBFROMUI(d64->bytes+4, targar[1]);
     }
    else {
-    pu[0]=targar[1];		   /* directly store the high int */
-    pu[1]=targar[0];		   /* then the low int */
+    /* hi int then lo */
+    UBFROMUI(d64->bytes,   targar[1]);
+    UBFROMUI(d64->bytes+4, targar[0]);
     }
 
   if (status!=0) decContextSetStatus(set, status); /* pass on status */
@@ -218,21 +219,20 @@ decNumber * decimal64ToNumber(const decimal64 *d64, decNumber *dn) {
   uInt msd;			   /* coefficient MSD */
   uInt exp;			   /* exponent top two bits */
   uInt comb;			   /* combination field */
-  const uInt *pu;		   /* work */
-  Int  need;			   /* .. */
+  Int  need;			   /* work */
+  uInt uiwork;			   /* for macros */
   uInt sourar[2];		   /* source 64-bit */
   #define sourhi sourar[1]	   /* name the word with the sign */
   #define sourlo sourar[0]	   /* and the lower word */
 
   /* load source from storage; this is endian */
-  pu=(const uInt *)d64->bytes;	   /* overlay */
   if (DECLITEND) {
-    sourlo=pu[0];		   /* directly load the low int */
-    sourhi=pu[1];		   /* then the high int */
+    sourlo=UBTOUI(d64->bytes  );   /* directly load the low int */
+    sourhi=UBTOUI(d64->bytes+4);   /* then the high int */
     }
    else {
-    sourhi=pu[0];		   /* directly load the high int */
-    sourlo=pu[1];		   /* then the low int */
+    sourhi=UBTOUI(d64->bytes  );   /* directly load the high int */
+    sourlo=UBTOUI(d64->bytes+4);   /* then the low int */
     }
 
   comb=(sourhi>>26)&0x1f;	   /* combination field */
@@ -307,25 +307,24 @@ char * decimal64ToString(const decimal64 *d64, char *string){
   uInt comb;			   /* combination field */
   char *cstart; 		   /* coefficient start */
   char *c;			   /* output pointer in string */
-  const uInt *pu;		   /* work */
+  const uByte *u;		   /* work */
   char *s, *t;			   /* .. (source, target) */
   Int  dpd;			   /* .. */
   Int  pre, e;			   /* .. */
-  const uByte *u;		   /* .. */
+  uInt uiwork;			   /* for macros */
 
   uInt sourar[2];		   /* source 64-bit */
   #define sourhi sourar[1]	   /* name the word with the sign */
   #define sourlo sourar[0]	   /* and the lower word */
 
   /* load source from storage; this is endian */
-  pu=(const uInt *)d64->bytes;	   /* overlay */
   if (DECLITEND) {
-    sourlo=pu[0];		   /* directly load the low int */
-    sourhi=pu[1];		   /* then the high int */
+    sourlo=UBTOUI(d64->bytes  );   /* directly load the low int */
+    sourhi=UBTOUI(d64->bytes+4);   /* then the high int */
     }
    else {
-    sourhi=pu[0];		   /* directly load the high int */
-    sourlo=pu[1];		   /* then the low int */
+    sourhi=UBTOUI(d64->bytes  );   /* directly load the high int */
+    sourlo=UBTOUI(d64->bytes+4);   /* then the low int */
     }
 
   c=string;			   /* where result will go */
@@ -472,7 +471,7 @@ decimal64 * decimal64FromString(decimal64 *result, const char *string,
 /*   returns 1 if the encoding of d64 is canonical, 0 otherwise       */
 /* No error is possible.					      */
 /* ------------------------------------------------------------------ */
-uint32_t decimal64IsCanonical(const decimal64 *d64) {
+uInt decimal64IsCanonical(const decimal64 *d64) {
   decNumber dn; 			/* work */
   decimal64 canon;			/* .. */
   decContext dc;			/* .. */
@@ -520,8 +519,8 @@ decimal64 * decimal64Canonical(decimal64 *result, const decimal64 *d64) {
 /* This assumes range has been checked and exponent previously 0; */
 /* type of exponent must be unsigned */
 #define decimal64SetExpCon(d, e) {				      \
-  (d)->bytes[0]|=(uint8_t)((e)>>6);				      \
-  (d)->bytes[1]|=(uint8_t)(((e)&0x3F)<<2);}
+  (d)->bytes[0]|=(uByte)((e)>>6);				      \
+  (d)->bytes[1]|=(uByte)(((e)&0x3F)<<2);}
 
 /* ------------------------------------------------------------------ */
 /* decimal64Show -- display a decimal64 in hexadecimal [debug aid]    */

libdecnumber/dpd/decimal64.h

@@ -82,7 +82,7 @@
   /* Routines							      */
   /* ---------------------------------------------------------------- */
 
-  #include "decimal64Symbols.h"
+#include "decimal64Symbols.h"
 
   /* String conversions 					      */
   decimal64 * decimal64FromString(decimal64 *, const char *, decContext *);