![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to fpu_log.c CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [Isaki's NoNo m68k/m88k emulator] / nono / fpe|
Return to fpu_log.c CVS log | Up to [Isaki's NoNo m68k/m88k emulator] / nono / fpe |
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.