[12c6f2d] | 1 | /*
|
---|
[df4ed85] | 2 | * Copyright (c) 2005 Josef Cejka
|
---|
[12c6f2d] | 3 | * All rights reserved.
|
---|
| 4 | *
|
---|
| 5 | * Redistribution and use in source and binary forms, with or without
|
---|
| 6 | * modification, are permitted provided that the following conditions
|
---|
| 7 | * are met:
|
---|
| 8 | *
|
---|
| 9 | * - Redistributions of source code must retain the above copyright
|
---|
| 10 | * notice, this list of conditions and the following disclaimer.
|
---|
| 11 | * - Redistributions in binary form must reproduce the above copyright
|
---|
| 12 | * notice, this list of conditions and the following disclaimer in the
|
---|
| 13 | * documentation and/or other materials provided with the distribution.
|
---|
| 14 | * - The name of the author may not be used to endorse or promote products
|
---|
| 15 | * derived from this software without specific prior written permission.
|
---|
| 16 | *
|
---|
| 17 | * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
|
---|
| 18 | * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
|
---|
| 19 | * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
|
---|
| 20 | * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
|
---|
| 21 | * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT
|
---|
| 22 | * NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
|
---|
| 23 | * DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
|
---|
| 24 | * THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
|
---|
| 25 | * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
|
---|
| 26 | * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
|
---|
| 27 | */
|
---|
| 28 |
|
---|
[231a60a] | 29 | /** @addtogroup softfloat
|
---|
[846848a6] | 30 | * @{
|
---|
| 31 | */
|
---|
| 32 | /** @file
|
---|
| 33 | */
|
---|
| 34 |
|
---|
[12c6f2d] | 35 | #include<sftypes.h>
|
---|
| 36 | #include<add.h>
|
---|
| 37 | #include<comparison.h>
|
---|
| 38 |
|
---|
| 39 | /** Add two Float32 numbers with same signs
|
---|
| 40 | */
|
---|
| 41 | float32 addFloat32(float32 a, float32 b)
|
---|
| 42 | {
|
---|
| 43 | int expdiff;
|
---|
[aa59fa0] | 44 | uint32_t exp1, exp2,frac1, frac2;
|
---|
[12c6f2d] | 45 |
|
---|
[4a5abddd] | 46 | expdiff = a.parts.exp - b.parts.exp;
|
---|
| 47 | if (expdiff < 0) {
|
---|
[12c6f2d] | 48 | if (isFloat32NaN(b)) {
|
---|
[1266543] | 49 | /* TODO: fix SigNaN */
|
---|
[12c6f2d] | 50 | if (isFloat32SigNaN(b)) {
|
---|
| 51 | };
|
---|
| 52 |
|
---|
| 53 | return b;
|
---|
| 54 | };
|
---|
| 55 |
|
---|
[4a5abddd] | 56 | if (b.parts.exp == FLOAT32_MAX_EXPONENT) {
|
---|
[12c6f2d] | 57 | return b;
|
---|
| 58 | }
|
---|
| 59 |
|
---|
[1266543] | 60 | frac1 = b.parts.fraction;
|
---|
[4a5abddd] | 61 | exp1 = b.parts.exp;
|
---|
[1266543] | 62 | frac2 = a.parts.fraction;
|
---|
[4a5abddd] | 63 | exp2 = a.parts.exp;
|
---|
| 64 | expdiff *= -1;
|
---|
[12c6f2d] | 65 | } else {
|
---|
[1266543] | 66 | if ((isFloat32NaN(a)) || (isFloat32NaN(b))) {
|
---|
| 67 | /* TODO: fix SigNaN */
|
---|
[4a5abddd] | 68 | if (isFloat32SigNaN(a) || isFloat32SigNaN(b)) {
|
---|
[12c6f2d] | 69 | };
|
---|
[1266543] | 70 | return (isFloat32NaN(a)?a:b);
|
---|
[12c6f2d] | 71 | };
|
---|
| 72 |
|
---|
[4a5abddd] | 73 | if (a.parts.exp == FLOAT32_MAX_EXPONENT) {
|
---|
[12c6f2d] | 74 | return a;
|
---|
| 75 | }
|
---|
| 76 |
|
---|
[1266543] | 77 | frac1 = a.parts.fraction;
|
---|
[4a5abddd] | 78 | exp1 = a.parts.exp;
|
---|
[1266543] | 79 | frac2 = b.parts.fraction;
|
---|
[4a5abddd] | 80 | exp2 = b.parts.exp;
|
---|
[12c6f2d] | 81 | };
|
---|
| 82 |
|
---|
| 83 | if (exp1 == 0) {
|
---|
| 84 | /* both are denormalized */
|
---|
[1266543] | 85 | frac1 += frac2;
|
---|
| 86 | if (frac1 & FLOAT32_HIDDEN_BIT_MASK ) {
|
---|
[12c6f2d] | 87 | /* result is not denormalized */
|
---|
| 88 | a.parts.exp = 1;
|
---|
| 89 | };
|
---|
[1266543] | 90 | a.parts.fraction = frac1;
|
---|
[12c6f2d] | 91 | return a;
|
---|
| 92 | };
|
---|
| 93 |
|
---|
[1266543] | 94 | frac1 |= FLOAT32_HIDDEN_BIT_MASK; /* add hidden bit */
|
---|
[12c6f2d] | 95 |
|
---|
| 96 | if (exp2 == 0) {
|
---|
| 97 | /* second operand is denormalized */
|
---|
[e6a40ac] | 98 | --expdiff;
|
---|
[12c6f2d] | 99 | } else {
|
---|
| 100 | /* add hidden bit to second operand */
|
---|
[1266543] | 101 | frac2 |= FLOAT32_HIDDEN_BIT_MASK;
|
---|
[12c6f2d] | 102 | };
|
---|
| 103 |
|
---|
| 104 | /* create some space for rounding */
|
---|
[1266543] | 105 | frac1 <<= 6;
|
---|
| 106 | frac2 <<= 6;
|
---|
[12c6f2d] | 107 |
|
---|
[1266543] | 108 | if (expdiff < (FLOAT32_FRACTION_SIZE + 2) ) {
|
---|
| 109 | frac2 >>= expdiff;
|
---|
| 110 | frac1 += frac2;
|
---|
[d3ca210] | 111 | } else {
|
---|
| 112 | a.parts.exp = exp1;
|
---|
| 113 | a.parts.fraction = (frac1 >> 6) & (~(FLOAT32_HIDDEN_BIT_MASK));
|
---|
| 114 | return a;
|
---|
| 115 | }
|
---|
[12c6f2d] | 116 |
|
---|
[1266543] | 117 | if (frac1 & (FLOAT32_HIDDEN_BIT_MASK << 7) ) {
|
---|
[12c6f2d] | 118 | ++exp1;
|
---|
[1266543] | 119 | frac1 >>= 1;
|
---|
[12c6f2d] | 120 | };
|
---|
| 121 |
|
---|
[1266543] | 122 | /* rounding - if first bit after fraction is set then round up */
|
---|
| 123 | frac1 += (0x1 << 5);
|
---|
[12c6f2d] | 124 |
|
---|
[1266543] | 125 | if (frac1 & (FLOAT32_HIDDEN_BIT_MASK << 7)) {
|
---|
| 126 | /* rounding overflow */
|
---|
[12c6f2d] | 127 | ++exp1;
|
---|
[1266543] | 128 | frac1 >>= 1;
|
---|
[12c6f2d] | 129 | };
|
---|
| 130 |
|
---|
[e6a40ac] | 131 |
|
---|
| 132 | if ((exp1 == FLOAT32_MAX_EXPONENT ) || (exp2 > exp1)) {
|
---|
[1266543] | 133 | /* overflow - set infinity as result */
|
---|
| 134 | a.parts.exp = FLOAT32_MAX_EXPONENT;
|
---|
| 135 | a.parts.fraction = 0;
|
---|
| 136 | return a;
|
---|
| 137 | }
|
---|
| 138 |
|
---|
[12c6f2d] | 139 | a.parts.exp = exp1;
|
---|
| 140 |
|
---|
| 141 | /*Clear hidden bit and shift */
|
---|
[1266543] | 142 | a.parts.fraction = ((frac1 >> 6) & (~FLOAT32_HIDDEN_BIT_MASK)) ;
|
---|
[12c6f2d] | 143 | return a;
|
---|
| 144 | }
|
---|
| 145 |
|
---|
| 146 |
|
---|
| 147 | /** Add two Float64 numbers with same signs
|
---|
| 148 | */
|
---|
| 149 | float64 addFloat64(float64 a, float64 b)
|
---|
| 150 | {
|
---|
| 151 | int expdiff;
|
---|
[aa59fa0] | 152 | uint32_t exp1, exp2;
|
---|
| 153 | uint64_t frac1, frac2;
|
---|
[12c6f2d] | 154 |
|
---|
[d3ca210] | 155 | expdiff = ((int )a.parts.exp) - b.parts.exp;
|
---|
[4a5abddd] | 156 | if (expdiff < 0) {
|
---|
[12c6f2d] | 157 | if (isFloat64NaN(b)) {
|
---|
[1266543] | 158 | /* TODO: fix SigNaN */
|
---|
[12c6f2d] | 159 | if (isFloat64SigNaN(b)) {
|
---|
| 160 | };
|
---|
| 161 |
|
---|
| 162 | return b;
|
---|
| 163 | };
|
---|
| 164 |
|
---|
| 165 | /* b is infinity and a not */
|
---|
[4a5abddd] | 166 | if (b.parts.exp == FLOAT64_MAX_EXPONENT ) {
|
---|
[12c6f2d] | 167 | return b;
|
---|
| 168 | }
|
---|
| 169 |
|
---|
[1266543] | 170 | frac1 = b.parts.fraction;
|
---|
[4a5abddd] | 171 | exp1 = b.parts.exp;
|
---|
[1266543] | 172 | frac2 = a.parts.fraction;
|
---|
[4a5abddd] | 173 | exp2 = a.parts.exp;
|
---|
| 174 | expdiff *= -1;
|
---|
[12c6f2d] | 175 | } else {
|
---|
| 176 | if (isFloat64NaN(a)) {
|
---|
[1266543] | 177 | /* TODO: fix SigNaN */
|
---|
[4a5abddd] | 178 | if (isFloat64SigNaN(a) || isFloat64SigNaN(b)) {
|
---|
[12c6f2d] | 179 | };
|
---|
| 180 | return a;
|
---|
| 181 | };
|
---|
| 182 |
|
---|
| 183 | /* a is infinity and b not */
|
---|
[4a5abddd] | 184 | if (a.parts.exp == FLOAT64_MAX_EXPONENT ) {
|
---|
[12c6f2d] | 185 | return a;
|
---|
| 186 | }
|
---|
| 187 |
|
---|
[1266543] | 188 | frac1 = a.parts.fraction;
|
---|
[12c6f2d] | 189 | exp1 = a.parts.exp;
|
---|
[1266543] | 190 | frac2 = b.parts.fraction;
|
---|
[12c6f2d] | 191 | exp2 = b.parts.exp;
|
---|
| 192 | };
|
---|
| 193 |
|
---|
| 194 | if (exp1 == 0) {
|
---|
| 195 | /* both are denormalized */
|
---|
[1266543] | 196 | frac1 += frac2;
|
---|
| 197 | if (frac1 & FLOAT64_HIDDEN_BIT_MASK) {
|
---|
[12c6f2d] | 198 | /* result is not denormalized */
|
---|
| 199 | a.parts.exp = 1;
|
---|
| 200 | };
|
---|
[1266543] | 201 | a.parts.fraction = frac1;
|
---|
[12c6f2d] | 202 | return a;
|
---|
| 203 | };
|
---|
| 204 |
|
---|
[1266543] | 205 | /* add hidden bit - frac1 is sure not denormalized */
|
---|
| 206 | frac1 |= FLOAT64_HIDDEN_BIT_MASK;
|
---|
[12c6f2d] | 207 |
|
---|
| 208 | /* second operand ... */
|
---|
| 209 | if (exp2 == 0) {
|
---|
| 210 | /* ... is denormalized */
|
---|
| 211 | --expdiff;
|
---|
| 212 | } else {
|
---|
| 213 | /* is not denormalized */
|
---|
[1266543] | 214 | frac2 |= FLOAT64_HIDDEN_BIT_MASK;
|
---|
[12c6f2d] | 215 | };
|
---|
| 216 |
|
---|
| 217 | /* create some space for rounding */
|
---|
[1266543] | 218 | frac1 <<= 6;
|
---|
| 219 | frac2 <<= 6;
|
---|
[12c6f2d] | 220 |
|
---|
[1266543] | 221 | if (expdiff < (FLOAT64_FRACTION_SIZE + 2) ) {
|
---|
| 222 | frac2 >>= expdiff;
|
---|
| 223 | frac1 += frac2;
|
---|
[d3ca210] | 224 | } else {
|
---|
| 225 | a.parts.exp = exp1;
|
---|
| 226 | a.parts.fraction = (frac1 >> 6) & (~(FLOAT64_HIDDEN_BIT_MASK));
|
---|
| 227 | return a;
|
---|
| 228 | }
|
---|
[12c6f2d] | 229 |
|
---|
[1266543] | 230 | if (frac1 & (FLOAT64_HIDDEN_BIT_MASK << 7) ) {
|
---|
[12c6f2d] | 231 | ++exp1;
|
---|
[1266543] | 232 | frac1 >>= 1;
|
---|
[12c6f2d] | 233 | };
|
---|
| 234 |
|
---|
[1266543] | 235 | /* rounding - if first bit after fraction is set then round up */
|
---|
| 236 | frac1 += (0x1 << 5);
|
---|
[12c6f2d] | 237 |
|
---|
[1266543] | 238 | if (frac1 & (FLOAT64_HIDDEN_BIT_MASK << 7)) {
|
---|
| 239 | /* rounding overflow */
|
---|
[12c6f2d] | 240 | ++exp1;
|
---|
[1266543] | 241 | frac1 >>= 1;
|
---|
[12c6f2d] | 242 | };
|
---|
| 243 |
|
---|
[d3ca210] | 244 | if ((exp1 == FLOAT64_MAX_EXPONENT ) || (exp2 > exp1)) {
|
---|
[1266543] | 245 | /* overflow - set infinity as result */
|
---|
| 246 | a.parts.exp = FLOAT64_MAX_EXPONENT;
|
---|
| 247 | a.parts.fraction = 0;
|
---|
| 248 | return a;
|
---|
| 249 | }
|
---|
| 250 |
|
---|
[12c6f2d] | 251 | a.parts.exp = exp1;
|
---|
| 252 | /*Clear hidden bit and shift */
|
---|
[1266543] | 253 | a.parts.fraction = ( (frac1 >> 6 ) & (~FLOAT64_HIDDEN_BIT_MASK));
|
---|
[d3ca210] | 254 |
|
---|
[12c6f2d] | 255 | return a;
|
---|
| 256 | }
|
---|
| 257 |
|
---|
[231a60a] | 258 | /** @}
|
---|
[846848a6] | 259 | */
|
---|