github.com/igggame/nebulas-go@v2.1.0+incompatible/nbre/3rd_party/SoftFloat-3e/source/s_mulAddF128.c (about) 1 2 /*============================================================================ 3 4 This C source file is part of the SoftFloat IEEE Floating-Point Arithmetic 5 Package, Release 3e, by John R. Hauser. 6 7 Copyright 2011, 2012, 2013, 2014, 2015, 2016 The Regents of the University of 8 California. All rights reserved. 9 10 Redistribution and use in source and binary forms, with or without 11 modification, are permitted provided that the following conditions are met: 12 13 1. Redistributions of source code must retain the above copyright notice, 14 this list of conditions, and the following disclaimer. 15 16 2. Redistributions in binary form must reproduce the above copyright notice, 17 this list of conditions, and the following disclaimer in the documentation 18 and/or other materials provided with the distribution. 19 20 3. Neither the name of the University nor the names of its contributors may 21 be used to endorse or promote products derived from this software without 22 specific prior written permission. 23 24 THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS "AS IS", AND ANY 25 EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED 26 WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE, ARE 27 DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE FOR ANY 28 DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES 29 (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; 30 LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND 31 ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 32 (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 33 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 34 35 =============================================================================*/ 36 37 #include <stdbool.h> 38 #include <stdint.h> 39 #include "platform.h" 40 #include "internals.h" 41 #include "specialize.h" 42 #include "softfloat.h" 43 44 float128_t 45 softfloat_mulAddF128( 46 uint_fast64_t uiA64, 47 uint_fast64_t uiA0, 48 uint_fast64_t uiB64, 49 uint_fast64_t uiB0, 50 uint_fast64_t uiC64, 51 uint_fast64_t uiC0, 52 uint_fast8_t op 53 ) 54 { 55 bool signA; 56 int_fast32_t expA; 57 struct uint128 sigA; 58 bool signB; 59 int_fast32_t expB; 60 struct uint128 sigB; 61 bool signC; 62 int_fast32_t expC; 63 struct uint128 sigC; 64 bool signZ; 65 uint_fast64_t magBits; 66 struct uint128 uiZ; 67 struct exp32_sig128 normExpSig; 68 int_fast32_t expZ; 69 uint64_t sig256Z[4]; 70 struct uint128 sigZ; 71 int_fast32_t shiftDist, expDiff; 72 struct uint128 x128; 73 uint64_t sig256C[4]; 74 static uint64_t zero256[4] = INIT_UINTM4( 0, 0, 0, 0 ); 75 uint_fast64_t sigZExtra, sig256Z0; 76 union ui128_f128 uZ; 77 78 /*------------------------------------------------------------------------ 79 *------------------------------------------------------------------------*/ 80 signA = signF128UI64( uiA64 ); 81 expA = expF128UI64( uiA64 ); 82 sigA.v64 = fracF128UI64( uiA64 ); 83 sigA.v0 = uiA0; 84 signB = signF128UI64( uiB64 ); 85 expB = expF128UI64( uiB64 ); 86 sigB.v64 = fracF128UI64( uiB64 ); 87 sigB.v0 = uiB0; 88 signC = signF128UI64( uiC64 ) ^ (op == softfloat_mulAdd_subC); 89 expC = expF128UI64( uiC64 ); 90 sigC.v64 = fracF128UI64( uiC64 ); 91 sigC.v0 = uiC0; 92 signZ = signA ^ signB ^ (op == softfloat_mulAdd_subProd); 93 /*------------------------------------------------------------------------ 94 *------------------------------------------------------------------------*/ 95 if ( expA == 0x7FFF ) { 96 if ( 97 (sigA.v64 | sigA.v0) || ((expB == 0x7FFF) && (sigB.v64 | sigB.v0)) 98 ) { 99 goto propagateNaN_ABC; 100 } 101 magBits = expB | sigB.v64 | sigB.v0; 102 goto infProdArg; 103 } 104 if ( expB == 0x7FFF ) { 105 if ( sigB.v64 | sigB.v0 ) goto propagateNaN_ABC; 106 magBits = expA | sigA.v64 | sigA.v0; 107 goto infProdArg; 108 } 109 if ( expC == 0x7FFF ) { 110 if ( sigC.v64 | sigC.v0 ) { 111 uiZ.v64 = 0; 112 uiZ.v0 = 0; 113 goto propagateNaN_ZC; 114 } 115 uiZ.v64 = uiC64; 116 uiZ.v0 = uiC0; 117 goto uiZ; 118 } 119 /*------------------------------------------------------------------------ 120 *------------------------------------------------------------------------*/ 121 if ( ! expA ) { 122 if ( ! (sigA.v64 | sigA.v0) ) goto zeroProd; 123 normExpSig = softfloat_normSubnormalF128Sig( sigA.v64, sigA.v0 ); 124 expA = normExpSig.exp; 125 sigA = normExpSig.sig; 126 } 127 if ( ! expB ) { 128 if ( ! (sigB.v64 | sigB.v0) ) goto zeroProd; 129 normExpSig = softfloat_normSubnormalF128Sig( sigB.v64, sigB.v0 ); 130 expB = normExpSig.exp; 131 sigB = normExpSig.sig; 132 } 133 /*------------------------------------------------------------------------ 134 *------------------------------------------------------------------------*/ 135 expZ = expA + expB - 0x3FFE; 136 sigA.v64 |= UINT64_C( 0x0001000000000000 ); 137 sigB.v64 |= UINT64_C( 0x0001000000000000 ); 138 sigA = softfloat_shortShiftLeft128( sigA.v64, sigA.v0, 8 ); 139 sigB = softfloat_shortShiftLeft128( sigB.v64, sigB.v0, 15 ); 140 softfloat_mul128To256M( sigA.v64, sigA.v0, sigB.v64, sigB.v0, sig256Z ); 141 sigZ.v64 = sig256Z[indexWord( 4, 3 )]; 142 sigZ.v0 = sig256Z[indexWord( 4, 2 )]; 143 shiftDist = 0; 144 if ( ! (sigZ.v64 & UINT64_C( 0x0100000000000000 )) ) { 145 --expZ; 146 shiftDist = -1; 147 } 148 if ( ! expC ) { 149 if ( ! (sigC.v64 | sigC.v0) ) { 150 shiftDist += 8; 151 goto sigZ; 152 } 153 normExpSig = softfloat_normSubnormalF128Sig( sigC.v64, sigC.v0 ); 154 expC = normExpSig.exp; 155 sigC = normExpSig.sig; 156 } 157 sigC.v64 |= UINT64_C( 0x0001000000000000 ); 158 sigC = softfloat_shortShiftLeft128( sigC.v64, sigC.v0, 8 ); 159 /*------------------------------------------------------------------------ 160 *------------------------------------------------------------------------*/ 161 expDiff = expZ - expC; 162 if ( expDiff < 0 ) { 163 expZ = expC; 164 if ( (signZ == signC) || (expDiff < -1) ) { 165 shiftDist -= expDiff; 166 if ( shiftDist ) { 167 sigZ = 168 softfloat_shiftRightJam128( sigZ.v64, sigZ.v0, shiftDist ); 169 } 170 } else { 171 if ( ! shiftDist ) { 172 x128 = 173 softfloat_shortShiftRight128( 174 sig256Z[indexWord( 4, 1 )], sig256Z[indexWord( 4, 0 )], 175 1 176 ); 177 sig256Z[indexWord( 4, 1 )] = (sigZ.v0<<63) | x128.v64; 178 sig256Z[indexWord( 4, 0 )] = x128.v0; 179 sigZ = softfloat_shortShiftRight128( sigZ.v64, sigZ.v0, 1 ); 180 sig256Z[indexWord( 4, 3 )] = sigZ.v64; 181 sig256Z[indexWord( 4, 2 )] = sigZ.v0; 182 } 183 } 184 } else { 185 if ( shiftDist ) softfloat_add256M( sig256Z, sig256Z, sig256Z ); 186 if ( ! expDiff ) { 187 sigZ.v64 = sig256Z[indexWord( 4, 3 )]; 188 sigZ.v0 = sig256Z[indexWord( 4, 2 )]; 189 } else { 190 sig256C[indexWord( 4, 3 )] = sigC.v64; 191 sig256C[indexWord( 4, 2 )] = sigC.v0; 192 sig256C[indexWord( 4, 1 )] = 0; 193 sig256C[indexWord( 4, 0 )] = 0; 194 softfloat_shiftRightJam256M( sig256C, expDiff, sig256C ); 195 } 196 } 197 /*------------------------------------------------------------------------ 198 *------------------------------------------------------------------------*/ 199 shiftDist = 8; 200 if ( signZ == signC ) { 201 /*-------------------------------------------------------------------- 202 *--------------------------------------------------------------------*/ 203 if ( expDiff <= 0 ) { 204 sigZ = softfloat_add128( sigC.v64, sigC.v0, sigZ.v64, sigZ.v0 ); 205 } else { 206 softfloat_add256M( sig256Z, sig256C, sig256Z ); 207 sigZ.v64 = sig256Z[indexWord( 4, 3 )]; 208 sigZ.v0 = sig256Z[indexWord( 4, 2 )]; 209 } 210 if ( sigZ.v64 & UINT64_C( 0x0200000000000000 ) ) { 211 ++expZ; 212 shiftDist = 9; 213 } 214 } else { 215 /*-------------------------------------------------------------------- 216 *--------------------------------------------------------------------*/ 217 if ( expDiff < 0 ) { 218 signZ = signC; 219 if ( expDiff < -1 ) { 220 sigZ = 221 softfloat_sub128( sigC.v64, sigC.v0, sigZ.v64, sigZ.v0 ); 222 sigZExtra = 223 sig256Z[indexWord( 4, 1 )] | sig256Z[indexWord( 4, 0 )]; 224 if ( sigZExtra ) { 225 sigZ = softfloat_sub128( sigZ.v64, sigZ.v0, 0, 1 ); 226 } 227 if ( ! (sigZ.v64 & UINT64_C( 0x0100000000000000 )) ) { 228 --expZ; 229 shiftDist = 7; 230 } 231 goto shiftRightRoundPack; 232 } else { 233 sig256C[indexWord( 4, 3 )] = sigC.v64; 234 sig256C[indexWord( 4, 2 )] = sigC.v0; 235 sig256C[indexWord( 4, 1 )] = 0; 236 sig256C[indexWord( 4, 0 )] = 0; 237 softfloat_sub256M( sig256C, sig256Z, sig256Z ); 238 } 239 } else if ( ! expDiff ) { 240 sigZ = softfloat_sub128( sigZ.v64, sigZ.v0, sigC.v64, sigC.v0 ); 241 if ( 242 ! (sigZ.v64 | sigZ.v0) && ! sig256Z[indexWord( 4, 1 )] 243 && ! sig256Z[indexWord( 4, 0 )] 244 ) { 245 goto completeCancellation; 246 } 247 sig256Z[indexWord( 4, 3 )] = sigZ.v64; 248 sig256Z[indexWord( 4, 2 )] = sigZ.v0; 249 if ( sigZ.v64 & UINT64_C( 0x8000000000000000 ) ) { 250 signZ = ! signZ; 251 softfloat_sub256M( zero256, sig256Z, sig256Z ); 252 } 253 } else { 254 softfloat_sub256M( sig256Z, sig256C, sig256Z ); 255 if ( 1 < expDiff ) { 256 sigZ.v64 = sig256Z[indexWord( 4, 3 )]; 257 sigZ.v0 = sig256Z[indexWord( 4, 2 )]; 258 if ( ! (sigZ.v64 & UINT64_C( 0x0100000000000000 )) ) { 259 --expZ; 260 shiftDist = 7; 261 } 262 goto sigZ; 263 } 264 } 265 /*-------------------------------------------------------------------- 266 *--------------------------------------------------------------------*/ 267 sigZ.v64 = sig256Z[indexWord( 4, 3 )]; 268 sigZ.v0 = sig256Z[indexWord( 4, 2 )]; 269 sigZExtra = sig256Z[indexWord( 4, 1 )]; 270 sig256Z0 = sig256Z[indexWord( 4, 0 )]; 271 if ( sigZ.v64 ) { 272 if ( sig256Z0 ) sigZExtra |= 1; 273 } else { 274 expZ -= 64; 275 sigZ.v64 = sigZ.v0; 276 sigZ.v0 = sigZExtra; 277 sigZExtra = sig256Z0; 278 if ( ! sigZ.v64 ) { 279 expZ -= 64; 280 sigZ.v64 = sigZ.v0; 281 sigZ.v0 = sigZExtra; 282 sigZExtra = 0; 283 if ( ! sigZ.v64 ) { 284 expZ -= 64; 285 sigZ.v64 = sigZ.v0; 286 sigZ.v0 = 0; 287 } 288 } 289 } 290 shiftDist = softfloat_countLeadingZeros64( sigZ.v64 ); 291 expZ += 7 - shiftDist; 292 shiftDist = 15 - shiftDist; 293 if ( 0 < shiftDist ) goto shiftRightRoundPack; 294 if ( shiftDist ) { 295 shiftDist = -shiftDist; 296 sigZ = softfloat_shortShiftLeft128( sigZ.v64, sigZ.v0, shiftDist ); 297 x128 = softfloat_shortShiftLeft128( 0, sigZExtra, shiftDist ); 298 sigZ.v0 |= x128.v64; 299 sigZExtra = x128.v0; 300 } 301 goto roundPack; 302 } 303 sigZ: 304 sigZExtra = sig256Z[indexWord( 4, 1 )] | sig256Z[indexWord( 4, 0 )]; 305 shiftRightRoundPack: 306 sigZExtra = (uint64_t) (sigZ.v0<<(64 - shiftDist)) | (sigZExtra != 0); 307 sigZ = softfloat_shortShiftRight128( sigZ.v64, sigZ.v0, shiftDist ); 308 roundPack: 309 return 310 softfloat_roundPackToF128( 311 signZ, expZ - 1, sigZ.v64, sigZ.v0, sigZExtra ); 312 /*------------------------------------------------------------------------ 313 *------------------------------------------------------------------------*/ 314 propagateNaN_ABC: 315 uiZ = softfloat_propagateNaNF128UI( uiA64, uiA0, uiB64, uiB0 ); 316 goto propagateNaN_ZC; 317 /*------------------------------------------------------------------------ 318 *------------------------------------------------------------------------*/ 319 infProdArg: 320 if ( magBits ) { 321 uiZ.v64 = packToF128UI64( signZ, 0x7FFF, 0 ); 322 uiZ.v0 = 0; 323 if ( expC != 0x7FFF ) goto uiZ; 324 if ( sigC.v64 | sigC.v0 ) goto propagateNaN_ZC; 325 if ( signZ == signC ) goto uiZ; 326 } 327 softfloat_raiseFlags( softfloat_flag_invalid ); 328 uiZ.v64 = defaultNaNF128UI64; 329 uiZ.v0 = defaultNaNF128UI0; 330 propagateNaN_ZC: 331 uiZ = softfloat_propagateNaNF128UI( uiZ.v64, uiZ.v0, uiC64, uiC0 ); 332 goto uiZ; 333 /*------------------------------------------------------------------------ 334 *------------------------------------------------------------------------*/ 335 zeroProd: 336 uiZ.v64 = uiC64; 337 uiZ.v0 = uiC0; 338 if ( ! (expC | sigC.v64 | sigC.v0) && (signZ != signC) ) { 339 completeCancellation: 340 uiZ.v64 = 341 packToF128UI64( 342 (softfloat_roundingMode == softfloat_round_min), 0, 0 ); 343 uiZ.v0 = 0; 344 } 345 uiZ: 346 uZ.ui = uiZ; 347 return uZ.f; 348 349 } 350