1/* IEEE754 floating point arithmetic 2 * double precision: common utilities 3 */ 4/* 5 * MIPS floating point support 6 * Copyright (C) 1994-2000 Algorithmics Ltd. 7 * http://www.algor.co.uk 8 * 9 * ######################################################################## 10 * 11 * This program is free software; you can distribute it and/or modify it 12 * under the terms of the GNU General Public License (Version 2) as 13 * published by the Free Software Foundation. 14 * 15 * This program is distributed in the hope it will be useful, but WITHOUT 16 * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 17 * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 18 * for more details. 19 * 20 * You should have received a copy of the GNU General Public License along 21 * with this program; if not, write to the Free Software Foundation, Inc., 22 * 59 Temple Place - Suite 330, Boston MA 02111-1307, USA. 23 * 24 * ######################################################################## 25 * 26 */ 27 28 29#include "ieee754dp.h" 30 31ieee754dp ieee754dp_add(ieee754dp x, ieee754dp y) 32{ 33 COMPXDP; 34 COMPYDP; 35 36 EXPLODEXDP; 37 EXPLODEYDP; 38 39 CLEARCX; 40 41 FLUSHXDP; 42 FLUSHYDP; 43 44 switch (CLPAIR(xc, yc)) { 45 case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_QNAN): 46 case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_SNAN): 47 case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_SNAN): 48 case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_SNAN): 49 case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_SNAN): 50 case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_SNAN): 51 case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_SNAN): 52 case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_ZERO): 53 case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_NORM): 54 case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_DNORM): 55 case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_INF): 56 SETCX(IEEE754_INVALID_OPERATION); 57 return ieee754dp_nanxcpt(ieee754dp_indef(), "add", x, y); 58 59 case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_QNAN): 60 case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_QNAN): 61 case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_QNAN): 62 case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_QNAN): 63 return y; 64 65 case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_QNAN): 66 case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_ZERO): 67 case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_NORM): 68 case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_DNORM): 69 case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_INF): 70 return x; 71 72 73 /* Infinity handling 74 */ 75 76 case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_INF): 77 if (xs == ys) 78 return x; 79 SETCX(IEEE754_INVALID_OPERATION); 80 return ieee754dp_xcpt(ieee754dp_indef(), "add", x, y); 81 82 case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_INF): 83 case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_INF): 84 case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_INF): 85 return y; 86 87 case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_ZERO): 88 case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_NORM): 89 case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_DNORM): 90 return x; 91 92 /* Zero handling 93 */ 94 95 case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_ZERO): 96 if (xs == ys) 97 return x; 98 else 99 return ieee754dp_zero(ieee754_csr.rm == 100 IEEE754_RD); 101 102 case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_ZERO): 103 case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_ZERO): 104 return x; 105 106 case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_NORM): 107 case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_DNORM): 108 return y; 109 110 case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_DNORM): 111 DPDNORMX; 112 113 /* FALL THROUGH */ 114 115 case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_DNORM): 116 DPDNORMY; 117 break; 118 119 case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_NORM): 120 DPDNORMX; 121 break; 122 123 case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_NORM): 124 break; 125 } 126 assert(xm & DP_HIDDEN_BIT); 127 assert(ym & DP_HIDDEN_BIT); 128 129 /* provide guard,round and stick bit space */ 130 xm <<= 3; 131 ym <<= 3; 132 133 if (xe > ye) { 134 /* have to shift y fraction right to align 135 */ 136 int s = xe - ye; 137 ym = XDPSRS(ym, s); 138 ye += s; 139 } else if (ye > xe) { 140 /* have to shift x fraction right to align 141 */ 142 int s = ye - xe; 143 xm = XDPSRS(xm, s); 144 xe += s; 145 } 146 assert(xe == ye); 147 assert(xe <= DP_EMAX); 148 149 if (xs == ys) { 150 /* generate 28 bit result of adding two 27 bit numbers 151 * leaving result in xm,xs,xe 152 */ 153 xm = xm + ym; 154 xe = xe; 155 xs = xs; 156 157 if (xm >> (DP_MBITS + 1 + 3)) { /* carry out */ 158 xm = XDPSRS1(xm); 159 xe++; 160 } 161 } else { 162 if (xm >= ym) { 163 xm = xm - ym; 164 xe = xe; 165 xs = xs; 166 } else { 167 xm = ym - xm; 168 xe = xe; 169 xs = ys; 170 } 171 if (xm == 0) 172 return ieee754dp_zero(ieee754_csr.rm == 173 IEEE754_RD); 174 175 /* normalize to rounding precision */ 176 while ((xm >> (DP_MBITS + 3)) == 0) { 177 xm <<= 1; 178 xe--; 179 } 180 181 } 182 DPNORMRET2(xs, xe, xm, "add", x, y); 183} 184