|
|
1.1 ! root 1: /* $NetBSD: fpu_log.c,v 1.18 2014/01/04 13:23:22 isaki Exp $ */ ! 2: ! 3: /* ! 4: * Copyright (c) 1995 Ken Nakata ! 5: * All rights reserved. ! 6: * ! 7: * Redistribution and use in source and binary forms, with or without ! 8: * modification, are permitted provided that the following conditions ! 9: * are met: ! 10: * 1. Redistributions of source code must retain the above copyright ! 11: * notice, this list of conditions and the following disclaimer. ! 12: * 2. Redistributions in binary form must reproduce the above copyright ! 13: * notice, this list of conditions and the following disclaimer in the ! 14: * documentation and/or other materials provided with the distribution. ! 15: * 3. Neither the name of the author nor the names of its contributors ! 16: * may be used to endorse or promote products derived from this software ! 17: * without specific prior written permission. ! 18: * ! 19: * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND ! 20: * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE ! 21: * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ! 22: * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE ! 23: * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL ! 24: * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS ! 25: * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) ! 26: * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT ! 27: * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY ! 28: * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF ! 29: * SUCH DAMAGE. ! 30: * ! 31: * @(#)fpu_log.c 10/8/95 ! 32: */ ! 33: ! 34: #include "fpu_emulate.h" ! 35: ! 36: static uint32_t logA6[] = { 0x3FC2499A, 0xB5E4040B }; ! 37: static uint32_t logA5[] = { 0xBFC555B5, 0x848CB7DB }; ! 38: static uint32_t logA4[] = { 0x3FC99999, 0x987D8730 }; ! 39: static uint32_t logA3[] = { 0xBFCFFFFF, 0xFF6F7E97 }; ! 40: static uint32_t logA2[] = { 0x3FD55555, 0x555555A4 }; ! 41: static uint32_t logA1[] = { 0xBFE00000, 0x00000008 }; ! 42: ! 43: static uint32_t logB5[] = { 0x3F175496, 0xADD7DAD6 }; ! 44: static uint32_t logB4[] = { 0x3F3C71C2, 0xFE80C7E0 }; ! 45: static uint32_t logB3[] = { 0x3F624924, 0x928BCCFF }; ! 46: static uint32_t logB2[] = { 0x3F899999, 0x999995EC }; ! 47: static uint32_t logB1[] = { 0x3FB55555, 0x55555555 }; ! 48: ! 49: /* sfpn = shortened fp number; can represent only positive numbers */ ! 50: static struct sfpn { ! 51: int sp_exp; ! 52: uint32_t sp_m0, sp_m1; ! 53: } logtbl[] = { ! 54: { 0x3FFE - 0x3fff, 0xFE03F80FU, 0xE03F80FEU }, ! 55: { 0x3FF7 - 0x3fff, 0xFF015358U, 0x833C47E2U }, ! 56: { 0x3FFE - 0x3fff, 0xFA232CF2U, 0x52138AC0U }, ! 57: { 0x3FF9 - 0x3fff, 0xBDC8D83EU, 0xAD88D549U }, ! 58: { 0x3FFE - 0x3fff, 0xF6603D98U, 0x0F6603DAU }, ! 59: { 0x3FFA - 0x3fff, 0x9CF43DCFU, 0xF5EAFD48U }, ! 60: { 0x3FFE - 0x3fff, 0xF2B9D648U, 0x0F2B9D65U }, ! 61: { 0x3FFA - 0x3fff, 0xDA16EB88U, 0xCB8DF614U }, ! 62: { 0x3FFE - 0x3fff, 0xEF2EB71FU, 0xC4345238U }, ! 63: { 0x3FFB - 0x3fff, 0x8B29B775U, 0x1BD70743U }, ! 64: { 0x3FFE - 0x3fff, 0xEBBDB2A5U, 0xC1619C8CU }, ! 65: { 0x3FFB - 0x3fff, 0xA8D839F8U, 0x30C1FB49U }, ! 66: { 0x3FFE - 0x3fff, 0xE865AC7BU, 0x7603A197U }, ! 67: { 0x3FFB - 0x3fff, 0xC61A2EB1U, 0x8CD907ADU }, ! 68: { 0x3FFE - 0x3fff, 0xE525982AU, 0xF70C880EU }, ! 69: { 0x3FFB - 0x3fff, 0xE2F2A47AU, 0xDE3A18AFU }, ! 70: { 0x3FFE - 0x3fff, 0xE1FC780EU, 0x1FC780E2U }, ! 71: { 0x3FFB - 0x3fff, 0xFF64898EU, 0xDF55D551U }, ! 72: { 0x3FFE - 0x3fff, 0xDEE95C4CU, 0xA037BA57U }, ! 73: { 0x3FFC - 0x3fff, 0x8DB956A9U, 0x7B3D0148U }, ! 74: { 0x3FFE - 0x3fff, 0xDBEB61EEU, 0xD19C5958U }, ! 75: { 0x3FFC - 0x3fff, 0x9B8FE100U, 0xF47BA1DEU }, ! 76: { 0x3FFE - 0x3fff, 0xD901B203U, 0x6406C80EU }, ! 77: { 0x3FFC - 0x3fff, 0xA9372F1DU, 0x0DA1BD17U }, ! 78: { 0x3FFE - 0x3fff, 0xD62B80D6U, 0x2B80D62CU }, ! 79: { 0x3FFC - 0x3fff, 0xB6B07F38U, 0xCE90E46BU }, ! 80: { 0x3FFE - 0x3fff, 0xD3680D36U, 0x80D3680DU }, ! 81: { 0x3FFC - 0x3fff, 0xC3FD0329U, 0x06488481U }, ! 82: { 0x3FFE - 0x3fff, 0xD0B69FCBU, 0xD2580D0BU }, ! 83: { 0x3FFC - 0x3fff, 0xD11DE0FFU, 0x15AB18CAU }, ! 84: { 0x3FFE - 0x3fff, 0xCE168A77U, 0x25080CE1U }, ! 85: { 0x3FFC - 0x3fff, 0xDE1433A1U, 0x6C66B150U }, ! 86: { 0x3FFE - 0x3fff, 0xCB8727C0U, 0x65C393E0U }, ! 87: { 0x3FFC - 0x3fff, 0xEAE10B5AU, 0x7DDC8ADDU }, ! 88: { 0x3FFE - 0x3fff, 0xC907DA4EU, 0x871146ADU }, ! 89: { 0x3FFC - 0x3fff, 0xF7856E5EU, 0xE2C9B291U }, ! 90: { 0x3FFE - 0x3fff, 0xC6980C69U, 0x80C6980CU }, ! 91: { 0x3FFD - 0x3fff, 0x82012CA5U, 0xA68206D7U }, ! 92: { 0x3FFE - 0x3fff, 0xC4372F85U, 0x5D824CA6U }, ! 93: { 0x3FFD - 0x3fff, 0x882C5FCDU, 0x7256A8C5U }, ! 94: { 0x3FFE - 0x3fff, 0xC1E4BBD5U, 0x95F6E947U }, ! 95: { 0x3FFD - 0x3fff, 0x8E44C60BU, 0x4CCFD7DEU }, ! 96: { 0x3FFE - 0x3fff, 0xBFA02FE8U, 0x0BFA02FFU }, ! 97: { 0x3FFD - 0x3fff, 0x944AD09EU, 0xF4351AF6U }, ! 98: { 0x3FFE - 0x3fff, 0xBD691047U, 0x07661AA3U }, ! 99: { 0x3FFD - 0x3fff, 0x9A3EECD4U, 0xC3EAA6B2U }, ! 100: { 0x3FFE - 0x3fff, 0xBB3EE721U, 0xA54D880CU }, ! 101: { 0x3FFD - 0x3fff, 0xA0218434U, 0x353F1DE8U }, ! 102: { 0x3FFE - 0x3fff, 0xB92143FAU, 0x36F5E02EU }, ! 103: { 0x3FFD - 0x3fff, 0xA5F2FCABU, 0xBBC506DAU }, ! 104: { 0x3FFE - 0x3fff, 0xB70FBB5AU, 0x19BE3659U }, ! 105: { 0x3FFD - 0x3fff, 0xABB3B8BAU, 0x2AD362A5U }, ! 106: { 0x3FFE - 0x3fff, 0xB509E68AU, 0x9B94821FU }, ! 107: { 0x3FFD - 0x3fff, 0xB1641795U, 0xCE3CA97BU }, ! 108: { 0x3FFE - 0x3fff, 0xB30F6352U, 0x8917C80BU }, ! 109: { 0x3FFD - 0x3fff, 0xB7047551U, 0x5D0F1C61U }, ! 110: { 0x3FFE - 0x3fff, 0xB11FD3B8U, 0x0B11FD3CU }, ! 111: { 0x3FFD - 0x3fff, 0xBC952AFEU, 0xEA3D13E1U }, ! 112: { 0x3FFE - 0x3fff, 0xAF3ADDC6U, 0x80AF3ADEU }, ! 113: { 0x3FFD - 0x3fff, 0xC2168ED0U, 0xF458BA4AU }, ! 114: { 0x3FFE - 0x3fff, 0xAD602B58U, 0x0AD602B6U }, ! 115: { 0x3FFD - 0x3fff, 0xC788F439U, 0xB3163BF1U }, ! 116: { 0x3FFE - 0x3fff, 0xAB8F69E2U, 0x8359CD11U }, ! 117: { 0x3FFD - 0x3fff, 0xCCECAC08U, 0xBF04565DU }, ! 118: { 0x3FFE - 0x3fff, 0xA9C84A47U, 0xA07F5638U }, ! 119: { 0x3FFD - 0x3fff, 0xD2420487U, 0x2DD85160U }, ! 120: { 0x3FFE - 0x3fff, 0xA80A80A8U, 0x0A80A80BU }, ! 121: { 0x3FFD - 0x3fff, 0xD7894992U, 0x3BC3588AU }, ! 122: { 0x3FFE - 0x3fff, 0xA655C439U, 0x2D7B73A8U }, ! 123: { 0x3FFD - 0x3fff, 0xDCC2C4B4U, 0x9887DACCU }, ! 124: { 0x3FFE - 0x3fff, 0xA4A9CF1DU, 0x96833751U }, ! 125: { 0x3FFD - 0x3fff, 0xE1EEBD3EU, 0x6D6A6B9EU }, ! 126: { 0x3FFE - 0x3fff, 0xA3065E3FU, 0xAE7CD0E0U }, ! 127: { 0x3FFD - 0x3fff, 0xE70D785CU, 0x2F9F5BDCU }, ! 128: { 0x3FFE - 0x3fff, 0xA16B312EU, 0xA8FC377DU }, ! 129: { 0x3FFD - 0x3fff, 0xEC1F392CU, 0x5179F283U }, ! 130: { 0x3FFE - 0x3fff, 0x9FD809FDU, 0x809FD80AU }, ! 131: { 0x3FFD - 0x3fff, 0xF12440D3U, 0xE36130E6U }, ! 132: { 0x3FFE - 0x3fff, 0x9E4CAD23U, 0xDD5F3A20U }, ! 133: { 0x3FFD - 0x3fff, 0xF61CCE92U, 0x346600BBU }, ! 134: { 0x3FFE - 0x3fff, 0x9CC8E160U, 0xC3FB19B9U }, ! 135: { 0x3FFD - 0x3fff, 0xFB091FD3U, 0x8145630AU }, ! 136: { 0x3FFE - 0x3fff, 0x9B4C6F9EU, 0xF03A3CAAU }, ! 137: { 0x3FFD - 0x3fff, 0xFFE97042U, 0xBFA4C2ADU }, ! 138: { 0x3FFE - 0x3fff, 0x99D722DAU, 0xBDE58F06U }, ! 139: { 0x3FFE - 0x3fff, 0x825EFCEDU, 0x49369330U }, ! 140: { 0x3FFE - 0x3fff, 0x9868C809U, 0x868C8098U }, ! 141: { 0x3FFE - 0x3fff, 0x84C37A7AU, 0xB9A905C9U }, ! 142: { 0x3FFE - 0x3fff, 0x97012E02U, 0x5C04B809U }, ! 143: { 0x3FFE - 0x3fff, 0x87224C2EU, 0x8E645FB7U }, ! 144: { 0x3FFE - 0x3fff, 0x95A02568U, 0x095A0257U }, ! 145: { 0x3FFE - 0x3fff, 0x897B8CACU, 0x9F7DE298U }, ! 146: { 0x3FFE - 0x3fff, 0x94458094U, 0x45809446U }, ! 147: { 0x3FFE - 0x3fff, 0x8BCF55DEU, 0xC4CD05FEU }, ! 148: { 0x3FFE - 0x3fff, 0x92F11384U, 0x0497889CU }, ! 149: { 0x3FFE - 0x3fff, 0x8E1DC0FBU, 0x89E125E5U }, ! 150: { 0x3FFE - 0x3fff, 0x91A2B3C4U, 0xD5E6F809U }, ! 151: { 0x3FFE - 0x3fff, 0x9066E68CU, 0x955B6C9BU }, ! 152: { 0x3FFE - 0x3fff, 0x905A3863U, 0x3E06C43BU }, ! 153: { 0x3FFE - 0x3fff, 0x92AADE74U, 0xC7BE59E0U }, ! 154: { 0x3FFE - 0x3fff, 0x8F1779D9U, 0xFDC3A219U }, ! 155: { 0x3FFE - 0x3fff, 0x94E9BFF6U, 0x15845643U }, ! 156: { 0x3FFE - 0x3fff, 0x8DDA5202U, 0x37694809U }, ! 157: { 0x3FFE - 0x3fff, 0x9723A1B7U, 0x20134203U }, ! 158: { 0x3FFE - 0x3fff, 0x8CA29C04U, 0x6514E023U }, ! 159: { 0x3FFE - 0x3fff, 0x995899C8U, 0x90EB8990U }, ! 160: { 0x3FFE - 0x3fff, 0x8B70344AU, 0x139BC75AU }, ! 161: { 0x3FFE - 0x3fff, 0x9B88BDAAU, 0x3A3DAE2FU }, ! 162: { 0x3FFE - 0x3fff, 0x8A42F870U, 0x5669DB46U }, ! 163: { 0x3FFE - 0x3fff, 0x9DB4224FU, 0xFFE1157CU }, ! 164: { 0x3FFE - 0x3fff, 0x891AC73AU, 0xE9819B50U }, ! 165: { 0x3FFE - 0x3fff, 0x9FDADC26U, 0x8B7A12DAU }, ! 166: { 0x3FFE - 0x3fff, 0x87F78087U, 0xF78087F8U }, ! 167: { 0x3FFE - 0x3fff, 0xA1FCFF17U, 0xCE733BD4U }, ! 168: { 0x3FFE - 0x3fff, 0x86D90544U, 0x7A34ACC6U }, ! 169: { 0x3FFE - 0x3fff, 0xA41A9E8FU, 0x5446FB9FU }, ! 170: { 0x3FFE - 0x3fff, 0x85BF3761U, 0x2CEE3C9BU }, ! 171: { 0x3FFE - 0x3fff, 0xA633CD7EU, 0x6771CD8BU }, ! 172: { 0x3FFE - 0x3fff, 0x84A9F9C8U, 0x084A9F9DU }, ! 173: { 0x3FFE - 0x3fff, 0xA8489E60U, 0x0B435A5EU }, ! 174: { 0x3FFE - 0x3fff, 0x83993052U, 0x3FBE3368U }, ! 175: { 0x3FFE - 0x3fff, 0xAA59233CU, 0xCCA4BD49U }, ! 176: { 0x3FFE - 0x3fff, 0x828CBFBEU, 0xB9A020A3U }, ! 177: { 0x3FFE - 0x3fff, 0xAC656DAEU, 0x6BCC4985U }, ! 178: { 0x3FFE - 0x3fff, 0x81848DA8U, 0xFAF0D277U }, ! 179: { 0x3FFE - 0x3fff, 0xAE6D8EE3U, 0x60BB2468U }, ! 180: { 0x3FFE - 0x3fff, 0x80808080U, 0x80808081U }, ! 181: { 0x3FFE - 0x3fff, 0xB07197A2U, 0x3C46C654U }, ! 182: }; ! 183: ! 184: static struct fpn *__fpu_logn(struct fpemu *fe); ! 185: ! 186: /* ! 187: * natural log - algorithm taken from Motorola FPSP, ! 188: * except this doesn't bother to check for invalid input. ! 189: */ ! 190: static struct fpn * ! 191: __fpu_logn(struct fpemu *fe) ! 192: { ! 193: static struct fpn X, F, U, V, W, KLOG2; ! 194: struct fpn *d; ! 195: int i, k; ! 196: ! 197: CPYFPN(&X, &fe->fe_f2); ! 198: ! 199: /* see if |X-1| < 1/16 approx. */ ! 200: if ((-1 == X.fp_exp && (0xf07d0000U >> (31 - FP_LG)) <= X.fp_mant[0]) || ! 201: (0 == X.fp_exp && X.fp_mant[0] <= (0x88410000U >> (31 - FP_LG)))) { ! 202: /* log near 1 */ ! 203: #if FPE_DEBUG ! 204: printf("__fpu_logn: log near 1\n"); ! 205: #endif ! 206: ! 207: fpu_const(&fe->fe_f1, FPU_CONST_1); ! 208: /* X+1 */ ! 209: d = fpu_add(fe); ! 210: CPYFPN(&V, d); ! 211: ! 212: CPYFPN(&fe->fe_f1, &X); ! 213: fpu_const(&fe->fe_f2, FPU_CONST_1); ! 214: fe->fe_f2.fp_sign = 1; /* -1.0 */ ! 215: /* X-1 */ ! 216: d = fpu_add(fe); ! 217: CPYFPN(&fe->fe_f1, d); ! 218: /* 2(X-1) */ ! 219: fe->fe_f1.fp_exp++; /* *= 2 */ ! 220: CPYFPN(&fe->fe_f2, &V); ! 221: /* U=2(X-1)/(X+1) */ ! 222: d = fpu_div(fe); ! 223: CPYFPN(&U, d); ! 224: CPYFPN(&fe->fe_f1, d); ! 225: CPYFPN(&fe->fe_f2, d); ! 226: /* V=U*U */ ! 227: d = fpu_mul(fe); ! 228: CPYFPN(&V, d); ! 229: CPYFPN(&fe->fe_f1, d); ! 230: CPYFPN(&fe->fe_f2, d); ! 231: /* W=V*V */ ! 232: d = fpu_mul(fe); ! 233: CPYFPN(&W, d); ! 234: ! 235: /* calculate U+U*V*([B1+W*(B3+W*B5)]+[V*(B2+W*B4)]) */ ! 236: ! 237: /* B1+W*(B3+W*B5) part */ ! 238: CPYFPN(&fe->fe_f1, d); ! 239: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logB5); ! 240: /* W*B5 */ ! 241: d = fpu_mul(fe); ! 242: CPYFPN(&fe->fe_f1, d); ! 243: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logB3); ! 244: /* B3+W*B5 */ ! 245: d = fpu_add(fe); ! 246: CPYFPN(&fe->fe_f1, d); ! 247: CPYFPN(&fe->fe_f2, &W); ! 248: /* W*(B3+W*B5) */ ! 249: d = fpu_mul(fe); ! 250: CPYFPN(&fe->fe_f1, d); ! 251: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logB1); ! 252: /* B1+W*(B3+W*B5) */ ! 253: d = fpu_add(fe); ! 254: CPYFPN(&X, d); ! 255: ! 256: /* [V*(B2+W*B4)] part */ ! 257: CPYFPN(&fe->fe_f1, &W); ! 258: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logB4); ! 259: /* W*B4 */ ! 260: d = fpu_mul(fe); ! 261: CPYFPN(&fe->fe_f1, d); ! 262: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logB2); ! 263: /* B2+W*B4 */ ! 264: d = fpu_add(fe); ! 265: CPYFPN(&fe->fe_f1, d); ! 266: CPYFPN(&fe->fe_f2, &V); ! 267: /* V*(B2+W*B4) */ ! 268: d = fpu_mul(fe); ! 269: CPYFPN(&fe->fe_f1, d); ! 270: CPYFPN(&fe->fe_f2, &X); ! 271: /* B1+W*(B3+W*B5)+V*(B2+W*B4) */ ! 272: d = fpu_add(fe); ! 273: CPYFPN(&fe->fe_f1, d); ! 274: CPYFPN(&fe->fe_f2, &V); ! 275: /* V*(B1+W*(B3+W*B5)+V*(B2+W*B4)) */ ! 276: d = fpu_mul(fe); ! 277: CPYFPN(&fe->fe_f1, d); ! 278: CPYFPN(&fe->fe_f2, &U); ! 279: /* U*V*(B1+W*(B3+W*B5)+V*(B2+W*B4)) */ ! 280: d = fpu_mul(fe); ! 281: CPYFPN(&fe->fe_f1, d); ! 282: CPYFPN(&fe->fe_f2, &U); ! 283: /* U+U*V*(B1+W*(B3+W*B5)+V*(B2+W*B4)) */ ! 284: d = fpu_add(fe); ! 285: } else /* the usual case */ { ! 286: #if FPE_DEBUG ! 287: printf("__fpu_logn: the usual case. X=(%d,%08x,%08x...)\n", ! 288: X.fp_exp, X.fp_mant[0], X.fp_mant[1]); ! 289: #endif ! 290: ! 291: k = X.fp_exp; ! 292: /* X <- Y */ ! 293: X.fp_exp = fe->fe_f2.fp_exp = 0; ! 294: ! 295: /* get the most significant 7 bits of X */ ! 296: F.fp_class = FPC_NUM; ! 297: F.fp_sign = 0; ! 298: F.fp_exp = X.fp_exp; ! 299: F.fp_mant[0] = X.fp_mant[0] & (0xfe000000U >> (31 - FP_LG)); ! 300: F.fp_mant[0] |= (0x01000000U >> (31 - FP_LG)); ! 301: F.fp_mant[1] = F.fp_mant[2] = 0; ! 302: F.fp_sticky = 0; ! 303: ! 304: #if FPE_DEBUG ! 305: printf("__fpu_logn: X=Y*2^k=(%d,%08x,%08x...)*2^%d\n", ! 306: fe->fe_f2.fp_exp, fe->fe_f2.fp_mant[0], ! 307: fe->fe_f2.fp_mant[1], k); ! 308: printf("__fpu_logn: F=(%d,%08x,%08x...)\n", ! 309: F.fp_exp, F.fp_mant[0], F.fp_mant[1]); ! 310: #endif ! 311: ! 312: /* index to the table */ ! 313: i = (F.fp_mant[0] >> (FP_LG - 7)) & 0x7e; ! 314: ! 315: #if FPE_DEBUG ! 316: printf("__fpu_logn: index to logtbl i=%d(%x)\n", i, i); ! 317: #endif ! 318: ! 319: CPYFPN(&fe->fe_f1, &F); ! 320: /* -F */ ! 321: fe->fe_f1.fp_sign = 1; ! 322: /* Y-F */ ! 323: d = fpu_add(fe); ! 324: CPYFPN(&fe->fe_f1, d); ! 325: ! 326: /* fe_f2 = 1/F */ ! 327: fe->fe_f2.fp_class = FPC_NUM; ! 328: fe->fe_f2.fp_sign = fe->fe_f2.fp_sticky = fe->fe_f2.fp_mant[2] ! 329: = 0; ! 330: fe->fe_f2.fp_exp = logtbl[i].sp_exp; ! 331: fe->fe_f2.fp_mant[0] = (logtbl[i].sp_m0 >> (31 - FP_LG)); ! 332: fe->fe_f2.fp_mant[1] = (logtbl[i].sp_m0 << (FP_LG + 1)) | ! 333: (logtbl[i].sp_m1 >> (31 - FP_LG)); ! 334: fe->fe_f2.fp_mant[2] = ! 335: (uint32_t)(logtbl[i].sp_m1 << (FP_LG + 1)); ! 336: ! 337: #if FPE_DEBUG ! 338: printf("__fpu_logn: 1/F=(%d,%08x,%08x...)\n", fe->fe_f2.fp_exp, ! 339: fe->fe_f2.fp_mant[0], fe->fe_f2.fp_mant[1]); ! 340: #endif ! 341: ! 342: /* U = (Y-F) * (1/F) */ ! 343: d = fpu_mul(fe); ! 344: CPYFPN(&U, d); ! 345: ! 346: /* KLOG2 = K * ln(2) */ ! 347: /* fe_f1 == (fpn)k */ ! 348: fpu_explode(fe, &fe->fe_f1, FTYPE_LNG, &k); ! 349: (void)fpu_const(&fe->fe_f2, FPU_CONST_LN_2); ! 350: #if FPE_DEBUG ! 351: printf("__fpu_logn: fp(k)=(%d,%08x,%08x...)\n", ! 352: fe->fe_f1.fp_exp, ! 353: fe->fe_f1.fp_mant[0], fe->fe_f1.fp_mant[1]); ! 354: printf("__fpu_logn: ln(2)=(%d,%08x,%08x...)\n", ! 355: fe->fe_f2.fp_exp, ! 356: fe->fe_f2.fp_mant[0], fe->fe_f2.fp_mant[1]); ! 357: #endif ! 358: /* K * LOGOF2 */ ! 359: d = fpu_mul(fe); ! 360: CPYFPN(&KLOG2, d); ! 361: ! 362: /* V=U*U */ ! 363: CPYFPN(&fe->fe_f1, &U); ! 364: CPYFPN(&fe->fe_f2, &U); ! 365: d = fpu_mul(fe); ! 366: CPYFPN(&V, d); ! 367: ! 368: /* ! 369: * approximation of LOG(1+U) by ! 370: * (U+V*(A1+V*(A3+V*A5)))+(U*V*(A2+V*(A4+V*A6))) ! 371: */ ! 372: ! 373: /* (U+V*(A1+V*(A3+V*A5))) part */ ! 374: CPYFPN(&fe->fe_f1, d); ! 375: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logA5); ! 376: /* V*A5 */ ! 377: d = fpu_mul(fe); ! 378: ! 379: CPYFPN(&fe->fe_f1, d); ! 380: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logA3); ! 381: /* A3+V*A5 */ ! 382: d = fpu_add(fe); ! 383: ! 384: CPYFPN(&fe->fe_f1, d); ! 385: CPYFPN(&fe->fe_f2, &V); ! 386: /* V*(A3+V*A5) */ ! 387: d = fpu_mul(fe); ! 388: ! 389: CPYFPN(&fe->fe_f1, d); ! 390: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logA1); ! 391: /* A1+V*(A3+V*A5) */ ! 392: d = fpu_add(fe); ! 393: ! 394: CPYFPN(&fe->fe_f1, d); ! 395: CPYFPN(&fe->fe_f2, &V); ! 396: /* V*(A1+V*(A3+V*A5)) */ ! 397: d = fpu_mul(fe); ! 398: ! 399: CPYFPN(&fe->fe_f1, d); ! 400: CPYFPN(&fe->fe_f2, &U); ! 401: /* U+V*(A1+V*(A3+V*A5)) */ ! 402: d = fpu_add(fe); ! 403: ! 404: CPYFPN(&X, d); ! 405: ! 406: /* (U*V*(A2+V*(A4+V*A6))) part */ ! 407: CPYFPN(&fe->fe_f1, &V); ! 408: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logA6); ! 409: /* V*A6 */ ! 410: d = fpu_mul(fe); ! 411: CPYFPN(&fe->fe_f1, d); ! 412: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logA4); ! 413: /* A4+V*A6 */ ! 414: d = fpu_add(fe); ! 415: CPYFPN(&fe->fe_f1, d); ! 416: CPYFPN(&fe->fe_f2, &V); ! 417: /* V*(A4+V*A6) */ ! 418: d = fpu_mul(fe); ! 419: CPYFPN(&fe->fe_f1, d); ! 420: fpu_explode(fe, &fe->fe_f2, FTYPE_DBL, logA2); ! 421: /* A2+V*(A4+V*A6) */ ! 422: d = fpu_add(fe); ! 423: CPYFPN(&fe->fe_f1, d); ! 424: CPYFPN(&fe->fe_f2, &V); ! 425: /* V*(A2+V*(A4+V*A6)) */ ! 426: d = fpu_mul(fe); ! 427: CPYFPN(&fe->fe_f1, d); ! 428: CPYFPN(&fe->fe_f2, &U); ! 429: /* U*V*(A2+V*(A4+V*A6)) */ ! 430: d = fpu_mul(fe); ! 431: CPYFPN(&fe->fe_f1, d); ! 432: i++; ! 433: /* fe_f2 = logtbl[i+1] (== LOG(F)) */ ! 434: fe->fe_f2.fp_class = FPC_NUM; ! 435: fe->fe_f2.fp_sign = fe->fe_f2.fp_sticky = fe->fe_f2.fp_mant[2] ! 436: = 0; ! 437: fe->fe_f2.fp_exp = logtbl[i].sp_exp; ! 438: fe->fe_f2.fp_mant[0] = (logtbl[i].sp_m0 >> (31 - FP_LG)); ! 439: fe->fe_f2.fp_mant[1] = (logtbl[i].sp_m0 << (FP_LG + 1)) | ! 440: (logtbl[i].sp_m1 >> (31 - FP_LG)); ! 441: fe->fe_f2.fp_mant[2] = (logtbl[i].sp_m1 << (FP_LG + 1)); ! 442: ! 443: #if FPE_DEBUG ! 444: printf("__fpu_logn: ln(F)=(%d,%08x,%08x,...)\n", ! 445: fe->fe_f2.fp_exp, ! 446: fe->fe_f2.fp_mant[0], fe->fe_f2.fp_mant[1]); ! 447: #endif ! 448: ! 449: /* LOG(F)+U*V*(A2+V*(A4+V*A6)) */ ! 450: d = fpu_add(fe); ! 451: CPYFPN(&fe->fe_f1, d); ! 452: CPYFPN(&fe->fe_f2, &X); ! 453: /* LOG(F)+U+V*(A1+V*(A3+V*A5))+U*V*(A2+V*(A4+V*A6)) */ ! 454: d = fpu_add(fe); ! 455: ! 456: #if FPE_DEBUG ! 457: printf("__fpu_logn: ln(Y)=(%c,%d,%08x,%08x,%08x)\n", ! 458: d->fp_sign ? '-' : '+', d->fp_exp, ! 459: d->fp_mant[0], d->fp_mant[1], d->fp_mant[2]); ! 460: #endif ! 461: ! 462: CPYFPN(&fe->fe_f1, d); ! 463: CPYFPN(&fe->fe_f2, &KLOG2); ! 464: /* K*LOGOF2+LOG(F)+U+V*(A1+V*(A3+V*A5))+U*V*(A2+V*(A4+V*A6)) */ ! 465: d = fpu_add(fe); ! 466: } ! 467: ! 468: return d; ! 469: } ! 470: ! 471: struct fpn * ! 472: fpu_log10(struct fpemu *fe) ! 473: { ! 474: struct fpn *fp = &fe->fe_f2; ! 475: uint32_t fpsr; ! 476: ! 477: fpsr = fe->fe_fpsr & ~FPSR_EXCP; /* clear all exceptions */ ! 478: ! 479: if (fp->fp_class >= FPC_NUM) { ! 480: if (fp->fp_sign) { /* negative number or Inf */ ! 481: fp = fpu_newnan(fe); ! 482: fpsr |= FPSR_OPERR; ! 483: } else if (fp->fp_class == FPC_NUM) { ! 484: /* the real work here */ ! 485: fp = __fpu_logn(fe); ! 486: if (fp != &fe->fe_f1) ! 487: CPYFPN(&fe->fe_f1, fp); ! 488: (void)fpu_const(&fe->fe_f2, FPU_CONST_LN_10); ! 489: fp = fpu_div(fe); ! 490: } /* else if fp == +Inf, return +Inf */ ! 491: } else if (fp->fp_class == FPC_ZERO) { ! 492: /* return -Inf */ ! 493: fp->fp_class = FPC_INF; ! 494: fp->fp_sign = 1; ! 495: fpsr |= FPSR_DZ; ! 496: } else if (fp->fp_class == FPC_SNAN) { ! 497: fpsr |= FPSR_SNAN; ! 498: fp = fpu_newnan(fe); ! 499: } else { ! 500: fp = fpu_newnan(fe); ! 501: } ! 502: ! 503: fe->fe_fpsr = fpsr; ! 504: ! 505: return fp; ! 506: } ! 507: ! 508: struct fpn * ! 509: fpu_log2(struct fpemu *fe) ! 510: { ! 511: struct fpn *fp = &fe->fe_f2; ! 512: uint32_t fpsr; ! 513: ! 514: fpsr = fe->fe_fpsr & ~FPSR_EXCP; /* clear all exceptions */ ! 515: ! 516: if (fp->fp_class >= FPC_NUM) { ! 517: if (fp->fp_sign) { /* negative number or Inf */ ! 518: fp = fpu_newnan(fe); ! 519: fpsr |= FPSR_OPERR; ! 520: } else if (fp->fp_class == FPC_NUM) { ! 521: /* the real work here */ ! 522: if (fp->fp_mant[0] == FP_1 && fp->fp_mant[1] == 0 && ! 523: fp->fp_mant[2] == 0) { ! 524: /* fp == 2.0 ^ exp <--> log2(fp) == exp */ ! 525: fpu_explode(fe, &fe->fe_f3, FTYPE_LNG, ! 526: &fp->fp_exp); ! 527: fp = &fe->fe_f3; ! 528: } else { ! 529: fp = __fpu_logn(fe); ! 530: if (fp != &fe->fe_f1) ! 531: CPYFPN(&fe->fe_f1, fp); ! 532: (void)fpu_const(&fe->fe_f2, FPU_CONST_LN_2); ! 533: fp = fpu_div(fe); ! 534: } ! 535: } /* else if fp == +Inf, return +Inf */ ! 536: } else if (fp->fp_class == FPC_ZERO) { ! 537: /* return -Inf */ ! 538: fp->fp_class = FPC_INF; ! 539: fp->fp_sign = 1; ! 540: fpsr |= FPSR_DZ; ! 541: } else if (fp->fp_class == FPC_SNAN) { ! 542: fpsr |= FPSR_SNAN; ! 543: fp = fpu_newnan(fe); ! 544: } else { ! 545: fp = fpu_newnan(fe); ! 546: } ! 547: ! 548: fe->fe_fpsr = fpsr; ! 549: return fp; ! 550: } ! 551: ! 552: struct fpn * ! 553: fpu_logn(struct fpemu *fe) ! 554: { ! 555: struct fpn *fp = &fe->fe_f2; ! 556: uint32_t fpsr; ! 557: ! 558: fpsr = fe->fe_fpsr & ~FPSR_EXCP; /* clear all exceptions */ ! 559: ! 560: if (fp->fp_class >= FPC_NUM) { ! 561: if (fp->fp_sign) { /* negative number or Inf */ ! 562: fp = fpu_newnan(fe); ! 563: fpsr |= FPSR_OPERR; ! 564: } else if (fp->fp_class == FPC_NUM) { ! 565: /* the real work here */ ! 566: fp = __fpu_logn(fe); ! 567: } /* else if fp == +Inf, return +Inf */ ! 568: } else if (fp->fp_class == FPC_ZERO) { ! 569: /* return -Inf */ ! 570: fp->fp_class = FPC_INF; ! 571: fp->fp_sign = 1; ! 572: fpsr |= FPSR_DZ; ! 573: } else if (fp->fp_class == FPC_SNAN) { ! 574: fpsr |= FPSR_SNAN; ! 575: fp = fpu_newnan(fe); ! 576: } else { ! 577: fp = fpu_newnan(fe); ! 578: } ! 579: ! 580: fe->fe_fpsr = fpsr; ! 581: ! 582: return fp; ! 583: } ! 584: ! 585: struct fpn * ! 586: fpu_lognp1(struct fpemu *fe) ! 587: { ! 588: struct fpn *fp; ! 589: ! 590: /* if src is +0/-0, return +0/-0 */ ! 591: if (ISZERO(&fe->fe_f2)) ! 592: return &fe->fe_f2; ! 593: ! 594: /* build a 1.0 */ ! 595: fp = fpu_const(&fe->fe_f1, FPU_CONST_1); ! 596: /* fp = 1.0 + f2 */ ! 597: fp = fpu_add(fe); ! 598: ! 599: /* copy the result to the src opr */ ! 600: CPYFPN(&fe->fe_f2, fp); ! 601: ! 602: return fpu_logn(fe); ! 603: }
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.