^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 1) // SPDX-License-Identifier: GPL-2.0-only
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 2) /* IEEE754 floating point arithmetic
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 3) * double precision: common utilities
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 4) */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 5) /*
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 6) * MIPS floating point support
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 7) * Copyright (C) 1994-2000 Algorithmics Ltd.
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 8) */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 9)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 10) #include "ieee754dp.h"
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 11)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 12) union ieee754dp ieee754dp_sub(union ieee754dp x, union ieee754dp y)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 13) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 14) int s;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 15)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 16) COMPXDP;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 17) COMPYDP;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 18)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 19) EXPLODEXDP;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 20) EXPLODEYDP;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 21)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 22) ieee754_clearcx();
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 23)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 24) FLUSHXDP;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 25) FLUSHYDP;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 26)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 27) switch (CLPAIR(xc, yc)) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 28) case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_SNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 29) case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_SNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 30) case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_SNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 31) case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_SNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 32) case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_SNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 33) return ieee754dp_nanxcpt(y);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 34)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 35) case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_SNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 36) case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_QNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 37) case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_ZERO):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 38) case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_NORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 39) case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_DNORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 40) case CLPAIR(IEEE754_CLASS_SNAN, IEEE754_CLASS_INF):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 41) return ieee754dp_nanxcpt(x);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 42)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 43) case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_QNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 44) case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_QNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 45) case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_QNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 46) case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_QNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 47) return y;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 48)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 49) case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_QNAN):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 50) case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_ZERO):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 51) case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_NORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 52) case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_DNORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 53) case CLPAIR(IEEE754_CLASS_QNAN, IEEE754_CLASS_INF):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 54) return x;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 55)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 56)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 57) /*
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 58) * Infinity handling
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 59) */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 60) case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_INF):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 61) if (xs != ys)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 62) return x;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 63) ieee754_setcx(IEEE754_INVALID_OPERATION);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 64) return ieee754dp_indef();
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 65)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 66) case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_INF):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 67) case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_INF):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 68) case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_INF):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 69) return ieee754dp_inf(ys ^ 1);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 70)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 71) case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_ZERO):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 72) case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_NORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 73) case CLPAIR(IEEE754_CLASS_INF, IEEE754_CLASS_DNORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 74) return x;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 75)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 76) /*
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 77) * Zero handling
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 78) */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 79) case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_ZERO):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 80) if (xs != ys)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 81) return x;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 82) else
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 83) return ieee754dp_zero(ieee754_csr.rm == FPU_CSR_RD);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 84)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 85) case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_ZERO):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 86) case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_ZERO):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 87) return x;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 88)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 89) case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_NORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 90) case CLPAIR(IEEE754_CLASS_ZERO, IEEE754_CLASS_DNORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 91) /* quick fix up */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 92) DPSIGN(y) ^= 1;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 93) return y;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 94)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 95) case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_DNORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 96) DPDNORMX;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 97) fallthrough;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 98) case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_DNORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 99) /* normalize ym,ye */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 100) DPDNORMY;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 101) break;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 102)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 103) case CLPAIR(IEEE754_CLASS_DNORM, IEEE754_CLASS_NORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 104) /* normalize xm,xe */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 105) DPDNORMX;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 106) break;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 107)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 108) case CLPAIR(IEEE754_CLASS_NORM, IEEE754_CLASS_NORM):
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 109) break;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 110) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 111) /* flip sign of y and handle as add */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 112) ys ^= 1;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 113)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 114) assert(xm & DP_HIDDEN_BIT);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 115) assert(ym & DP_HIDDEN_BIT);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 116)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 117)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 118) /* provide guard,round and stick bit dpace */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 119) xm <<= 3;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 120) ym <<= 3;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 121)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 122) if (xe > ye) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 123) /*
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 124) * Have to shift y fraction right to align
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 125) */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 126) s = xe - ye;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 127) ym = XDPSRS(ym, s);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 128) ye += s;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 129) } else if (ye > xe) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 130) /*
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 131) * Have to shift x fraction right to align
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 132) */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 133) s = ye - xe;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 134) xm = XDPSRS(xm, s);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 135) xe += s;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 136) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 137) assert(xe == ye);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 138) assert(xe <= DP_EMAX);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 139)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 140) if (xs == ys) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 141) /* generate 28 bit result of adding two 27 bit numbers
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 142) */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 143) xm = xm + ym;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 144)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 145) if (xm >> (DP_FBITS + 1 + 3)) { /* carry out */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 146) xm = XDPSRS1(xm); /* shift preserving sticky */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 147) xe++;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 148) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 149) } else {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 150) if (xm >= ym) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 151) xm = xm - ym;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 152) } else {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 153) xm = ym - xm;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 154) xs = ys;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 155) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 156) if (xm == 0) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 157) if (ieee754_csr.rm == FPU_CSR_RD)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 158) return ieee754dp_zero(1); /* round negative inf. => sign = -1 */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 159) else
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 160) return ieee754dp_zero(0); /* other round modes => sign = 1 */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 161) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 162)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 163) /* normalize to rounding precision
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 164) */
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 165) while ((xm >> (DP_FBITS + 3)) == 0) {
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 166) xm <<= 1;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 167) xe--;
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 168) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 169) }
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 170)
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 171) return ieee754dp_format(xs, xe, xm);
^8f3ce5b39 (kx 2023-10-28 12:00:06 +0300 172) }