[BACK]Return to fpu_log.c CVS log [TXT][DIR] Up to [local] / sys / arch / m68k / fpe

Annotation of sys/arch/m68k/fpe/fpu_log.c, Revision 1.1.1.1

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

CVSweb