usr/src/lib/libbc/libc/gen/common/_Qfaddsub.c

   1 /*
   2  * CDDL HEADER START
   3  *
   4  * The contents of this file are subject to the terms of the
   5  * Common Development and Distribution License, Version 1.0 only
   6  * (the "License").  You may not use this file except in compliance
   7  * with the License.
   8  *
   9  * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
  10  * or http://www.opensolaris.org/os/licensing.
  11  * See the License for the specific language governing permissions
  12  * and limitations under the License.
  13  *
  14  * When distributing Covered Code, include this CDDL HEADER in each
  15  * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
  16  * If applicable, add the following below this CDDL HEADER, with the
  17  * fields enclosed by brackets "[]" replaced with your own identifying
  18  * information: Portions Copyright [yyyy] [name of copyright owner]
  19  *
  20  * CDDL HEADER END
  21  */
  22 #pragma ident   "%Z%%M% %I%     %E% SMI"
  23
  24 /*
  25  * Copyright (c) 1988 by Sun Microsystems, Inc.
  26  */
  27
  28 #include "_Qquad.h"
  29 #include "_Qglobals.h"
  30
  31 PRIVATE void
  32 true_add(px, py, pz)
  33         unpacked       *px, *py, *pz;
  34
  35 {
  36         unsigned        c;
  37         unpacked       *pt;
  38
  39         if ((int) px->fpclass < (int) py->fpclass) {    /* Reverse. */
  40                 pt = py;
  41                 py = px;
  42                 px = pt;
  43         }
  44         /* Now class(x) >= class(y). */
  45         switch (px->fpclass) {
  46         case fp_quiet:          /* NaN + x -> NaN */
  47         case fp_signaling:      /* NaN + x -> NaN */
  48         case fp_infinity:       /* Inf + x -> Inf */
  49         case fp_zero:           /* 0 + 0 -> 0 */
  50                 *pz = *px;
  51                 return;
  52         default:
  53                 if (py->fpclass == fp_zero) {
  54                         *pz = *px;
  55                         return;
  56                 }
  57         }
  58         /* Now z is normal or subnormal. */
  59         /* Now y is normal or subnormal. */
  60         if (px->exponent < py->exponent) {      /* Reverse. */
  61                 pt = py;
  62                 py = px;
  63                 px = pt;
  64         }
  65         /* Now class(x) >= class(y). */
  66         pz->fpclass = px->fpclass;
  67         pz->sign = px->sign;
  68         pz->exponent = px->exponent;
  69         pz->rounded = pz->sticky  = 0;
  70
  71         if (px->exponent != py->exponent) {     /* pre-alignment required */
  72                 fpu_rightshift(py, pz->exponent - py->exponent);
  73                 pz->rounded = py->rounded;
  74                 pz->sticky  = py->sticky;
  75         }
  76         c = 0;
  77         c = fpu_add3wc(&(pz->significand[3]),px->significand[3],
  78                                                 py->significand[3],c);
  79         c = fpu_add3wc(&(pz->significand[2]),px->significand[2],
  80                                                 py->significand[2],c);
  81         c = fpu_add3wc(&(pz->significand[1]),px->significand[1],
  82                                                 py->significand[1],c);
  83         c = fpu_add3wc(&(pz->significand[0]),px->significand[0],
  84                                                 py->significand[0],c);
  85
  86         /* Handle carry out of msb. */
  87         if(pz->significand[0]>=0x20000) {
  88                 fpu_rightshift(pz, 1);  /* Carried out bit. */
  89                 pz->exponent ++;        /* Renormalize. */
  90         }
  91         return;
  92 }
  93
  94 PRIVATE void
  95 true_sub(px, py, pz)
  96         unpacked       *px, *py, *pz;
  97
  98 {
  99         unsigned       *z,g,s,r,c;
 100         int            n;
 101         unpacked       *pt;
 102
 103         if ((int) px->fpclass < (int) py->fpclass) {    /* Reverse. */
 104                 pt = py;
 105                 py = px;
 106                 px = pt;
 107         }
 108         /* Now class(x) >= class(y). */
 109         *pz = *px;              /* Tentative difference: x. */
 110         switch (pz->fpclass) {
 111         case fp_quiet:          /* NaN - x -> NaN */
 112         case fp_signaling:      /* NaN - x -> NaN */
 113                 return;
 114         case fp_infinity:       /* Inf - x -> Inf */
 115                 if (py->fpclass == fp_infinity) {
 116                         fpu_error_nan(pz);      /* Inf - Inf -> NaN */
 117                         pz->fpclass = fp_quiet;
 118                 }
 119                 return;
 120         case fp_zero:           /* 0 - 0 -> 0 */
 121                 pz->sign = (fp_direction == fp_negative);
 122                 return;
 123         default:
 124                 if (py->fpclass == fp_zero)
 125                         return;
 126         }
 127
 128         /* x and y are both normal or subnormal. */
 129
 130         if (px->exponent < py->exponent) { /* Reverse. */
 131                 pt = py;
 132                 py = px;
 133                 px = pt;
 134         }
 135         /* Now exp(x) >= exp(y). */
 136         pz->fpclass = px->fpclass;
 137         pz->sign = px->sign;
 138         pz->exponent = px->exponent;
 139         pz->rounded = 0;
 140         pz->sticky = 0;
 141         z = pz->significand;
 142
 143         if (px->exponent == py->exponent) {     /* no pre-alignment required */
 144                 c = 0;
 145                 c = fpu_sub3wc(&z[3],px->significand[3],py->significand[3],c);
 146                 c = fpu_sub3wc(&z[2],px->significand[2],py->significand[2],c);
 147                 c = fpu_sub3wc(&z[1],px->significand[1],py->significand[1],c);
 148                 c = fpu_sub3wc(&z[0],px->significand[0],py->significand[0],c);
 149                 if((z[0]|z[1]|z[2]|z[3])==0) {          /* exact zero result */
 150                         pz->sign = (fp_direction == fp_negative);
 151                         pz->fpclass = fp_zero;
 152                         return;
 153                 }
 154                 if(z[0]>=0x20000) {     /* sign reversal occurred */
 155                         pz->sign = py->sign;
 156                         c = 0;
 157                         c = fpu_neg2wc(&z[3],z[3],c);
 158                         c = fpu_neg2wc(&z[2],z[2],c);
 159                         c = fpu_neg2wc(&z[1],z[1],c);
 160                         c = fpu_neg2wc(&z[0],z[0],c);
 161                 }
 162                 fpu_normalize(pz);
 163                 return;
 164         } else {                /* pre-alignment required */
 165                 fpu_rightshift(py, pz->exponent - py->exponent - 1);
 166                 r = py->rounded;        /* rounded bit */
 167                 s = py->sticky;         /* sticky bit */
 168                 fpu_rightshift(py, 1);
 169                 g = py->rounded;        /* guard bit */
 170                 if(s!=0) r = (r==0);
 171                 if((r|s)!=0) g = (g==0);/* guard and rounded bits of z */
 172                 c = ((g|r|s)!=0);
 173                 c = fpu_sub3wc(&z[3],px->significand[3],py->significand[3],c);
 174                 c = fpu_sub3wc(&z[2],px->significand[2],py->significand[2],c);
 175                 c = fpu_sub3wc(&z[1],px->significand[1],py->significand[1],c);
 176                 c = fpu_sub3wc(&z[0],px->significand[0],py->significand[0],c);
 177
 178                 if(z[0]>=0x10000) {     /* don't need post-shifted */
 179                         pz->sticky = s|r;
 180                         pz->rounded = g;
 181                 } else {                /* post-shifted left 1 bit */
 182                         pz->sticky = s;
 183                         pz->rounded = r;
 184                         pz->significand[0] = (z[0]<<1)|((z[1]&0x80000000)>>31);
 185                         pz->significand[1] = (z[1]<<1)|((z[2]&0x80000000)>>31);
 186                         pz->significand[2] = (z[2]<<1)|((z[3]&0x80000000)>>31);
 187                         pz->significand[3] = (z[3]<<1)|g;
 188                         pz->exponent      -= 1;
 189                         if(z[0]<0x10000) fpu_normalize(pz);
 190                 }
 191                 return;
 192         }
 193 }
 194
 195 void
 196 _fp_add(px, py, pz)
 197         unpacked       *px, *py, *pz;
 198
 199 {
 200         if (px->sign == py->sign)
 201                 true_add(px, py, pz);
 202         else
 203                 true_sub(px, py, pz);
 204 }
 205
 206 void
 207 _fp_sub(px, py, pz)
 208         unpacked       *px, *py, *pz;
 209
 210 {
 211         py->sign = 1 - py->sign;
 212         if (px->sign == py->sign)
 213                 true_add(px, py, pz);
 214         else
 215                 true_sub(px, py, pz);
 216 }