![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to fpu_implode.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_implode.c CVS log | Up to [Isaki's NoNo m68k/m88k emulator] / nono / fpe |
1.1 root 1: /* $NetBSD: fpu_implode.c,v 1.15 2013/03/26 11:30:21 isaki Exp $ */
2:
3: /*
4: * Copyright (c) 1992, 1993
5: * The Regents of the University of California. All rights reserved.
6: *
7: * This software was developed by the Computer Systems Engineering group
8: * at Lawrence Berkeley Laboratory under DARPA contract BG 91-66 and
9: * contributed to Berkeley.
10: *
11: * All advertising materials mentioning features or use of this software
12: * must display the following acknowledgement:
13: * This product includes software developed by the University of
14: * California, Lawrence Berkeley Laboratory.
15: *
16: * Redistribution and use in source and binary forms, with or without
17: * modification, are permitted provided that the following conditions
18: * are met:
19: * 1. Redistributions of source code must retain the above copyright
20: * notice, this list of conditions and the following disclaimer.
21: * 2. Redistributions in binary form must reproduce the above copyright
22: * notice, this list of conditions and the following disclaimer in the
23: * documentation and/or other materials provided with the distribution.
24: * 3. Neither the name of the University nor the names of its contributors
25: * may be used to endorse or promote products derived from this software
26: * without specific prior written permission.
27: *
28: * THIS SOFTWARE IS PROVIDED BY THE REGENTS AND CONTRIBUTORS ``AS IS'' AND
29: * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
30: * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
31: * ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE
32: * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
33: * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
34: * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
35: * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
36: * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
37: * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
38: * SUCH DAMAGE.
39: *
40: * @(#)fpu_implode.c 8.1 (Berkeley) 6/11/93
41: */
42:
43: /*
44: * FPU subroutines: `implode' internal format numbers into the machine's
45: * `packed binary' format.
46: */
47:
48: #include "fpu_emulate.h"
49:
50: /* Conversion from internal format -- note asymmetry. */
51: static uint32_t fpu_ftoi(struct fpemu *fe, struct fpn *fp);
52: static uint32_t fpu_ftos(struct fpemu *fe, struct fpn *fp);
53: static uint32_t fpu_ftod(struct fpemu *fe, struct fpn *fp, uint32_t *);
54: static uint32_t fpu_ftox(struct fpemu *fe, struct fpn *fp, uint32_t *);
55: static void fpu_round_chkmin(struct fpemu *fe, struct fpn *fp, int, int);
56: static void fpu_round_chkinf(struct fpemu *fe, struct fpn *fp, int);
57:
58: /*
59: * Round a number (algorithm from Motorola MC68882 manual, modified for
60: * our internal format). Set inexact exception if rounding is required.
61: * Return true iff we rounded up.
62: *
63: * After rounding, we discard the guard and round bits by shifting right
64: * 2 bits (a la fpu_shr(), but we do not bother with fp->fp_sticky).
65: * This saves effort later.
66: *
67: * Note that we may leave the value 2.0 in fp->fp_mant; it is the caller's
68: * responsibility to fix this if necessary.
69: */
70: // FPCR_ROUND (NR, RZ, RP, RM) によって fp (拡張精度) をラウンディングする。
71: // 戻り値はラウンドアップしたら 1(true)。
72: // ラウンディングが発生すると内部 FPSR:EXCP に INEX2 を立てる。
73: // 任意精度のラウンディングは fpu_round_prec() 参照のこと。XXX どうしたもんか
74: //
75: // 実行後の fp_mant[2] は GR ビットを捨ててあり、つまり fp_mant 全体が当初
76: // に比べて2ビット右シフトされた状態となる。
77: // ラウンドアップによって生じた整数部 2 は何もせずそのまま返されるので、
78: // これを必要に応じてどうにかするのは呼び出し側の責任。
79: // (ただし整数部(小数点位置)も2ビットずれてるはず、でいいのかな?)
80: // fp_sticky は常にゼロクリアされる。
81: // fp_exp は変更しない。
82: int
83: fpu_round(struct fpemu *fe, struct fpn *fp)
84: {
85: uint32_t m0, m1, m2;
86: int gr, s;
87:
88: m0 = fp->fp_mant[0];
89: m1 = fp->fp_mant[1];
90: m2 = fp->fp_mant[2];
91: gr = m2 & 3;
92: s = fp->fp_sticky;
93:
94: /* mant >>= FP_NG */
95: m2 = (m2 >> FP_NG) | (m1 << (32 - FP_NG));
96: m1 = (m1 >> FP_NG) | (m0 << (32 - FP_NG));
97: m0 >>= FP_NG;
98:
99: if ((gr | s) == 0) /* result is exact: no rounding needed */
100: goto rounddown;
101:
102: fe->fe_fpsr |= FPSR_INEX2; /* inexact */
103:
104: /* Go to rounddown to round down; break to round up. */
105: switch (fe->fe_fpcr & FPCR_ROUND) {
106:
107: case FPCR_NEAR:
108: default:
109: /*
110: * Round only if guard is set (gr & 2). If guard is set,
111: * but round & sticky both clear, then we want to round
112: * but have a tie, so round to even, i.e., add 1 iff odd.
113: */
114: if ((gr & 2) == 0)
115: goto rounddown;
116: if ((gr & 1) || fp->fp_sticky || (m2 & 1))
117: break;
118: goto rounddown;
119:
120: case FPCR_ZERO:
121: /* Round towards zero, i.e., down. */
122: goto rounddown;
123:
124: case FPCR_MINF:
125: /* Round towards -Inf: up if negative, down if positive. */
126: if (fp->fp_sign)
127: break;
128: goto rounddown;
129:
130: case FPCR_PINF:
131: /* Round towards +Inf: up if positive, down otherwise. */
132: if (!fp->fp_sign)
133: break;
134: goto rounddown;
135: }
136:
137: /* Bump low bit of mantissa, with carry. */
138: if (++m2 == 0 && ++m1 == 0)
139: m0++;
140: fp->fp_sticky = 0;
141: fp->fp_mant[0] = m0;
142: fp->fp_mant[1] = m1;
143: fp->fp_mant[2] = m2;
144: return (1);
145:
146: rounddown:
147: fp->fp_sticky = 0;
148: fp->fp_mant[0] = m0;
149: fp->fp_mant[1] = m1;
150: fp->fp_mant[2] = m2;
151: return (0);
152: }
153:
154: /*
155: * For overflow: return true if overflow is to go to +/-Inf, according
156: * to the sign of the overflowing result. If false, overflow is to go
157: * to the largest magnitude value instead.
158: */
159: static int
160: toinf(struct fpemu *fe, int sign)
161: {
162: int inf;
163:
164: /* look at rounding direction */
165: switch (fe->fe_fpcr & FPCR_ROUND) {
166:
167: default:
168: case FPCR_NEAR: /* the nearest value is always Inf */
169: inf = 1;
170: break;
171:
172: case FPCR_ZERO: /* toward 0 => never towards Inf */
173: inf = 0;
174: break;
175:
176: case FPCR_PINF: /* toward +Inf iff positive */
177: inf = (sign == 0);
178: break;
179:
180: case FPCR_MINF: /* toward -Inf iff negative */
181: inf = sign;
182: break;
183: }
184: return (inf);
185: }
186:
187: /*
188: * fpn -> int (int value returned as return value).
189: *
190: * N.B.: this conversion always rounds towards zero (this is a peculiarity
191: * of the SPARC instruction set).
192: */
193: static uint32_t
194: fpu_ftoi(struct fpemu *fe, struct fpn *fp)
195: {
196: uint32_t i;
197: int sign, exp;
198:
199: sign = fp->fp_sign;
200: switch (fp->fp_class) {
201: case FPC_ZERO:
202: return (0);
203:
204: case FPC_NUM:
205: /*
206: * If exp >= 2^32, overflow. Otherwise shift value right
207: * into last mantissa word (this will not exceed 0xffffffff),
208: * shifting any guard and round bits out into the sticky
209: * bit. Then ``round'' towards zero, i.e., just set an
210: * inexact exception if sticky is set (see fpu_round()).
211: * If the result is > 0x80000000, or is positive and equals
212: * 0x80000000, overflow; otherwise the last fraction word
213: * is the result.
214: */
215: if ((exp = fp->fp_exp) >= 32)
216: break;
217: /* NB: the following includes exp < 0 cases */
218: if (fpu_shr(fp, FP_NMANT - 1 - FP_NG - exp) != 0) {
219: /*
220: * m68881/2 do not underflow when
221: * converting to integer
222: */
223: ;
224: }
225: fpu_round(fe, fp);
226: i = fp->fp_mant[2];
227: if (i >= ((uint32_t)0x80000000 + sign))
228: break;
229: return (sign ? -i : i);
230:
231: default: /* Inf, qNaN, sNaN */
232: break;
233: }
234: /* overflow: replace any inexact exception with invalid */
235: fe->fe_fpsr = (fe->fe_fpsr & ~FPSR_INEX2) | FPSR_OPERR;
236: return (0x7fffffff + sign);
237: }
238:
239: /*
240: * fpn -> single (32 bit single returned as return value).
241: * We assume <= 29 bits in a single-precision fraction (1.f part).
242: */
243: static uint32_t
244: fpu_ftos(struct fpemu *fe, struct fpn *fp)
245: {
246: uint32_t sign = fp->fp_sign << 31;
247: int exp;
248:
249: #define SNG_EXP(e) ((e) << SNG_FRACBITS) /* makes e an exponent */
250: #define SNG_MASK (SNG_EXP(1) - 1) /* mask for fraction */
251:
252: /* Take care of non-numbers first. */
253: if (ISNAN(fp)) {
254: /*
255: * Preserve upper bits of NaN, per SPARC V8 appendix N.
256: * Note that fp->fp_mant[0] has the quiet bit set,
257: * even if it is classified as a signalling NaN.
258: */
259: (void) fpu_shr(fp, FP_NMANT - 1 - SNG_FRACBITS);
260: exp = SNG_EXP_INFNAN;
261: goto done;
262: }
263: if (ISINF(fp))
264: return (sign | SNG_EXP(SNG_EXP_INFNAN));
265: if (ISZERO(fp))
266: return (sign);
267:
268: /*
269: * Normals (including subnormals). Drop all the fraction bits
270: * (including the explicit ``implied'' 1 bit) down into the
271: * single-precision range. If the number is subnormal, move
272: * the ``implied'' 1 into the explicit range as well, and shift
273: * right to introduce leading zeroes. Rounding then acts
274: * differently for normals and subnormals: the largest subnormal
275: * may round to the smallest normal (1.0 x 2^minexp), or may
276: * remain subnormal. In the latter case, signal an underflow
277: * if the result was inexact or if underflow traps are enabled.
278: *
279: * Rounding a normal, on the other hand, always produces another
280: * normal (although either way the result might be too big for
281: * single precision, and cause an overflow). If rounding a
282: * normal produces 2.0 in the fraction, we need not adjust that
283: * fraction at all, since both 1.0 and 2.0 are zero under the
284: * fraction mask.
285: *
286: * Note that the guard and round bits vanish from the number after
287: * rounding.
288: */
289: if ((exp = fp->fp_exp + SNG_EXP_BIAS) <= 0) { /* subnormal */
290: fe->fe_fpsr |= FPSR_UNFL;
291: /* -NG for g,r; -SNG_FRACBITS-exp for fraction */
292: (void) fpu_shr(fp, FP_NMANT - FP_NG - SNG_FRACBITS - exp);
293: if (fpu_round(fe, fp) && fp->fp_mant[2] == SNG_EXP(1))
294: return (sign | SNG_EXP(1) | 0);
295: if (fe->fe_fpsr & FPSR_INEX2) {
296: /* mc68881/2 don't underflow when converting */
297: fe->fe_fpsr |= FPSR_UNFL;
298: }
299: return (sign | SNG_EXP(0) | fp->fp_mant[2]);
300: }
301: /* -FP_NG for g,r; -1 for implied 1; -SNG_FRACBITS for fraction */
302: (void) fpu_shr(fp, FP_NMANT - FP_NG - 1 - SNG_FRACBITS);
303: #ifdef DIAGNOSTIC
304: if ((fp->fp_mant[2] & SNG_EXP(1 << FP_NG)) == 0)
305: panic("fpu_ftos");
306: #endif
307: if (fpu_round(fe, fp) && fp->fp_mant[2] == SNG_EXP(2))
308: exp++;
309: if (exp >= SNG_EXP_INFNAN) {
310: /* overflow to inf or to max single */
311: fe->fe_fpsr |= FPSR_OVFL;
312: if (toinf(fe, sign)) {
313: fp->fp_class = FPC_INF;
314: return (sign | SNG_EXP(SNG_EXP_INFNAN));
315: }
316: fe->fe_fpsr |= FPSR_OPERR;
317: return (sign | SNG_EXP(SNG_EXP_INFNAN - 1) | SNG_MASK);
318: }
319: done:
320: if ((fp->fp_mant[2] & SNG_MASK) == 0)
321: fp->fp_class = FPC_ZERO;
322: /* phew, made it */
323: return (sign | SNG_EXP(exp) | (fp->fp_mant[2] & SNG_MASK));
324: }
325:
326: /*
327: * fpn -> double (32 bit high-order result returned; 32-bit low order result
328: * left in res[1]). Assumes <= 61 bits in double precision fraction.
329: *
330: * This code mimics fpu_ftos; see it for comments.
331: */
332: static uint32_t
333: fpu_ftod(struct fpemu *fe, struct fpn *fp, uint32_t *res)
334: {
335: uint32_t sign = fp->fp_sign << 31;
336: int exp;
337:
338: #define DBL_EXP(e) ((e) << (DBL_FRACBITS & 31))
339: #define DBL_MASK (DBL_EXP(1) - 1)
340:
341: if (ISNAN(fp)) {
342: (void) fpu_shr(fp, FP_NMANT - 1 - DBL_FRACBITS);
343: exp = DBL_EXP_INFNAN;
344: goto done;
345: }
346: if (ISINF(fp)) {
347: sign |= DBL_EXP(DBL_EXP_INFNAN);
348: res[1] = 0;
349: return (sign);
350: }
351: if (ISZERO(fp)) {
352: res[1] = 0;
353: return (sign);
354: }
355:
356: if ((exp = fp->fp_exp + DBL_EXP_BIAS) <= 0) {
357: fe->fe_fpsr |= FPSR_UNFL;
358: (void) fpu_shr(fp, FP_NMANT - FP_NG - DBL_FRACBITS - exp);
359: if (fpu_round(fe, fp) && fp->fp_mant[1] == DBL_EXP(1)) {
360: res[1] = 0;
361: return (sign | DBL_EXP(1) | 0);
362: }
363: if (fe->fe_fpsr & FPSR_INEX2) {
364: /* mc68881/2 don't underflow when converting */
365: fe->fe_fpsr |= FPSR_UNFL;
366: }
367: exp = 0;
368: goto done;
369: }
370: (void) fpu_shr(fp, FP_NMANT - FP_NG - 1 - DBL_FRACBITS);
371: if (fpu_round(fe, fp) && fp->fp_mant[1] == DBL_EXP(2))
372: exp++;
373: if (exp >= DBL_EXP_INFNAN) {
374: fe->fe_fpsr |= FPSR_OVFL;
375: if (toinf(fe, sign)) {
376: fp->fp_class = FPC_INF;
377: res[1] = 0;
378: return (sign | DBL_EXP(DBL_EXP_INFNAN) | 0);
379: }
380: fe->fe_fpsr |= FPSR_OPERR;
381: res[1] = ~0;
382: return (sign | DBL_EXP(DBL_EXP_INFNAN) | DBL_MASK);
383: }
384: done:
385: res[1] = fp->fp_mant[2];
386: if (((fp->fp_mant[1] & DBL_MASK) | res[1]) == 0)
387: fp->fp_class = FPC_ZERO;
388: return (sign | DBL_EXP(exp) | (fp->fp_mant[1] & DBL_MASK));
389: }
390:
391: /*
392: * fpn -> 68k extended (32 bit high-order result returned; two 32-bit low
393: * order result left in res[1] & res[2]). Assumes == 64 bits in extended
394: * precision fraction.
395: *
396: * This code mimics fpu_ftos; see it for comments.
397: */
398: static uint32_t
399: fpu_ftox(struct fpemu *fe, struct fpn *fp, uint32_t *res)
400: {
401: uint32_t sign = fp->fp_sign << 31;
402: int exp;
403:
404: #define EXT_EXP(e) ((e) << 16)
405: /*
406: * on m68k extended prec, significand does not share the same long
407: * word with exponent
408: */
409: #define EXT_MASK 0
410: #define EXT_EXPLICIT1 (1UL << (63 & 31))
411: #define EXT_EXPLICIT2 (1UL << (64 & 31))
412:
413: if (ISNAN(fp)) {
414: (void) fpu_shr(fp, FP_NMANT - EXT_FRACBITS);
415: exp = EXT_EXP_INFNAN;
416: goto done;
417: }
418: if (ISINF(fp)) {
419: sign |= EXT_EXP(EXT_EXP_INFNAN);
420: res[1] = res[2] = 0;
421: return (sign);
422: }
423: if (ISZERO(fp)) {
424: res[1] = res[2] = 0;
425: return (sign);
426: }
427:
428: if ((exp = fp->fp_exp + EXT_EXP_BIAS) < 0) {
429: fe->fe_fpsr |= FPSR_UNFL;
430: /*
431: * I'm not sure about this <=... exp==0 doesn't mean
432: * it's a denormal in extended format
433: */
434: (void) fpu_shr(fp, FP_NMANT - FP_NG - EXT_FRACBITS - exp);
435: if (fpu_round(fe, fp) && fp->fp_mant[1] == EXT_EXPLICIT1) {
436: res[1] = res[2] = 0;
437: return (sign | EXT_EXP(1) | 0);
438: }
439: if (fe->fe_fpsr & FPSR_INEX2) {
440: /* mc68881/2 don't underflow */
441: fe->fe_fpsr |= FPSR_UNFL;
442: }
443: exp = 0;
444: goto done;
445: }
446: #if (FP_NMANT - FP_NG - EXT_FRACBITS) > 0
447: (void) fpu_shr(fp, FP_NMANT - FP_NG - EXT_FRACBITS);
448: #endif
449: if (fpu_round(fe, fp) && fp->fp_mant[0] == EXT_EXPLICIT2) {
450: exp++;
451: fpu_shr(fp, 1);
452: }
453: if (exp >= EXT_EXP_INFNAN) {
454: fe->fe_fpsr |= FPSR_OVFL;
455: if (toinf(fe, sign)) {
456: fp->fp_class = FPC_INF;
457: res[1] = res[2] = 0;
458: return (sign | EXT_EXP(EXT_EXP_INFNAN) | 0);
459: }
460: fe->fe_fpsr |= FPSR_OPERR;
461: res[1] = res[2] = ~0;
462: return (sign | EXT_EXP(EXT_EXP_INFNAN) | EXT_MASK);
463: }
464: done:
465: res[1] = fp->fp_mant[1];
466: res[2] = fp->fp_mant[2];
467: if ((res[1] | res[2]) == 0)
468: fp->fp_class = FPC_ZERO;
469: return (sign | EXT_EXP(exp));
470: }
471:
472: /*
473: * Implode an fpn, writing the result into the given space.
474: */
475: void
476: fpu_implode(struct fpemu *fe, struct fpn *fp, int type, uint32_t *space)
477: {
478: /* XXX Dont delete exceptions set here: fe->fe_fpsr &= ~FPSR_EXCP; */
479:
480: switch (type) {
481: case FTYPE_LNG:
482: space[0] = fpu_ftoi(fe, fp);
483: break;
484:
485: case FTYPE_SNG:
486: space[0] = fpu_ftos(fe, fp);
487: break;
488:
489: case FTYPE_DBL:
490: space[0] = fpu_ftod(fe, fp, space);
491: break;
492:
493: case FTYPE_EXT:
494: /* funky rounding precision options ?? */
495: space[0] = fpu_ftox(fe, fp, space);
496: break;
497:
498: default:
499: /* 何も出来ることがない */
500: break;
501: }
502: }
503:
504: #if defined(XM6i_FPE)
505: /*
506: * Shift the given number left lsh bits.
507: * Note that the sticky filed is cleared.
508: */
509: static void
510: fpu_shl(struct fpn *fp, int lsh)
511: {
512: uint32_t m0, m1, m2;
513: int rsh;
514:
515: m0 = fp->fp_mant[0];
516: m1 = fp->fp_mant[1];
517: m2 = fp->fp_mant[2];
518:
519: while (lsh >= 32) {
520: m0 = m1;
521: m1 = m2;
522: m2 = (fp->fp_sticky != 0) ? 0x80000000 : 0;
523: fp->fp_sticky = 0;
524: lsh -= 32;
525: }
526: if (lsh != 0) {
527: rsh = 32 - lsh;
528: m0 = (m0 << lsh) | (m1 >> rsh);
529: m1 = (m1 << lsh) | (m2 >> rsh);
530: m2 = (m2 << lsh);
531: m2 |= (fp->fp_sticky != 0) ? (1 << (lsh - 1)) : 0;
532: fp->fp_sticky = 0;
533: }
534: fp->fp_mant[0] = m0 & (FP_2 - 1);
535: fp->fp_mant[1] = m1;
536: fp->fp_mant[2] = m2;
537: }
538:
539: /*
540: * Round a number according to FPCR_MODE.
541: * (fpu_round() rounds according to FPCR_ROUND as FPCR_MODE = FPCR_EXTD)
542: */
543: // fp を FPCR_PREC(精度)/ FPCR_MODE(RN,RZ,RP,RM) によってラウンディングする。
544: // FPSR:EXCP に UNFL、INEX2 を立てる場合がある。
545: // 元からある fpu_round() は拡張精度限定のラウンディング。どうしたもんか。
546: void
547: fpu_round_prec(struct fpemu *fe, struct fpn *fp)
548: {
549: if (fp->fp_class != FPC_NUM) {
550: PRINTF("fpu_round_prec class != NUM\n");
551: return;
552: }
553:
554: switch ((fe->fe_fpcr & FPCR_PREC)) {
555: case FPCR_SNGL:
556: fpu_round_chkmin(fe, fp, SNG_EXP_BIAS, SNG_FRACBITS);
557: fpu_round_chkinf(fe, fp, SNG_EXP_BIAS);
558: break;
559:
560: case FPCR_DBL:
561: fpu_round_chkmin(fe, fp, DBL_EXP_BIAS, DBL_FRACBITS);
562: fpu_round_chkinf(fe, fp, DBL_EXP_BIAS);
563: break;
564:
565: case FPCR_EXTD:
566: default:
567: //
568: // 拡張精度だけ fpu_round_chk{min,inf}() と微妙に処理が異なる…。
569: //
570: if (fp->fp_exp < -EXT_EXP_BIAS - EXT_FRACBITS) {
571: // 拡張精度にすると指数が小さすぎて表現できない場合
572: fe->fe_fpsr |= FPSR_UNFL;
573: DUMPFP("e1:start", fp);
574:
575: // sticky だけの状態から FPCR_MODE によるラウンディングで
576: // fp_mant だけ作る。ここで仮数部ゼロならゼロ。
577: fp->fp_mant[0] = 0;
578: fp->fp_mant[1] = 0;
579: fp->fp_mant[2] = 0;
580: fp->fp_sticky = 1;
581: fpu_round(fe, fp);
582: if (fp->fp_mant[2] == 0)
583: fp->fp_class = FPC_ZERO;
584: // fp_mant の LSB に立ってるはずのビットを整数部になるよう正規化。
585: fp->fp_exp = FP_NMANT - EXT_EXP_BIAS - EXT_FRACBITS;
586: fpu_norm(fp);
587: } else if (fp->fp_exp < -EXT_EXP_BIAS) {
588: // 拡張精度にすると非正規化数になる場合
589: fe->fe_fpsr |= FPSR_UNFL;
590: DUMPFP("e2:start", fp);
591:
592: int shift = FP_NMANT - EXT_FRACBITS - 0/*integer bit included*/;
593: int effbits = EXT_EXP_BIAS + EXT_FRACBITS + fp->fp_exp;
594: PRINTF("effbits=%d\n", effbits);
595: shift += EXT_FRACBITS - effbits;
596: PRINTF("shift=%d\n", shift);
597: fpu_shr(fp, shift - FP_NG);
598: DUMPFP("e2:shr ", fp);
599:
600: int rup = fpu_round(fe, fp);
601: uint64_t m = (((uint64_t)fp->fp_mant[1]) << 32)
602: | (uint64_t)fp->fp_mant[2];
603: if (rup && m == (1ULL << effbits)) {
604: fp->fp_exp++;
605: shift--;
606: }
607: DUMPFP("e2:round", fp);
608:
609: fpu_shl(fp, shift);
610: DUMPFP("e2:shl ", fp);
611:
612: if (fp->fp_exp <= -EXT_EXP_BIAS - EXT_FRACBITS) {
613: PRINTF("e2:zero\n");
614: fp->fp_class = FPC_ZERO;
615: }
616: } else {
617: int shift = FP_NMANT - EXT_FRACBITS - 0/*Integer bit included*/;
618: PRINTF("shift=%d\n", shift);
619: fpu_shr(fp, shift - FP_NG);
620: DUMPFP("e3:shr ", fp);
621:
622: int rup = fpu_round(fe, fp);
623: uint32_t m = fp->fp_mant[1] | fp->fp_mant[2];
624: if (rup && fp->fp_mant[0] == 1 && m == 0) { /* mant == 2.0 */
625: fp->fp_exp++;
626: shift--;
627: DUMPFP("e3:rup ", fp);
628: } else
629: DUMPFP("e3:round", fp);
630:
631: fpu_shl(fp, shift);
632: DUMPFP("e3:shl ", fp);
633: }
634:
635: if (fp->fp_exp > EXT_EXP_BIAS) {
636: PRINTF("fpu_round_prec EXTD\n");
637: fe->fe_fpsr |= FPSR_OVFL;
638:
639: // 拡張精度で表現できない大きい値の場合は
640: // 正で RN/RP か、負で RN/RM なら Inf、
641: // そうでなければ拡張精度で表現できる最大値。
642: if (fp->fp_sign == 0) {
643: switch ((fe->fe_fpcr & FPCR_ROUND)) {
644: case FPCR_NEAR:
645: case FPCR_PINF:
646: fp->fp_class = FPC_INF;
647: break;
648: }
649: } else {
650: switch ((fe->fe_fpcr & FPCR_ROUND)) {
651: case FPCR_NEAR:
652: case FPCR_MINF:
653: fp->fp_class = FPC_INF;
654: break;
655: }
656: }
657: if (fp->fp_class != FPC_INF) {
658: PRINTF("fpu_round_prec make INF\n");
659: // 最大値
660: fp->fp_exp = EXT_EXP_BIAS;
661: fp->fp_mant[0] = FP_2 -1;
662: fp->fp_mant[1] = 0xffffffff;
663: fp->fp_mant[2] = 0xfff80000;
664: }
665: }
666: break;
667: }
668: }
669:
670: void
671: fpu_round_chkmin(struct fpemu *fe, struct fpn *fp, int EXP_BIAS, int FRACBITS)
672: {
673: int shift;
674: int effbits;
675: int rup;
676:
677: if (fp->fp_exp < -EXP_BIAS - FRACBITS) {
678: // 指定精度にすると指数が小さすぎて表現できない場合
679: fe->fe_fpsr |= FPSR_UNFL;
680: DUMPFP("s1:start", fp);
681:
682: // sticky だけの状態から FPCR_MODE によるラウンディングで
683: // fp_mant だけ作る。ここで仮数部ゼロならゼロ。
684: fp->fp_mant[0] = 0;
685: fp->fp_mant[1] = 0;
686: fp->fp_mant[2] = 0;
687: fp->fp_sticky = 1;
688: fpu_round(fe, fp);
689: DUMPFP("s1:round", fp);
690: if (fp->fp_mant[2] == 0)
691: fp->fp_class = FPC_ZERO;
692: // fp_mant の LSB に立ってるはずのビットを整数部になるよう正規化。
693: fp->fp_exp = FP_NMANT - EXP_BIAS - FRACBITS;
694: fpu_norm(fp);
695: DUMPFP("s1:norm ", fp);
696: return;
697: }
698: if (fp->fp_exp <= -EXP_BIAS) {
699: // 指定精度にすると非正規化数になる場合
700: fe->fe_fpsr |= FPSR_UNFL;
701: DUMPFP("s2:start", fp);
702:
703: shift = FP_NMANT - FRACBITS - 0/*No integer bit*/;
704: /*
705: * 1.XX~XX * 2^-127 => effbits = 22
706: * 1.X * 2^-148 => effbits = 2
707: * 1. * 2^-149 => effbits = 1
708: */
709: effbits = EXP_BIAS + FRACBITS + fp->fp_exp;
710: PRINTF("effbits=%d\n", effbits);
711: shift += FRACBITS - effbits;
712: PRINTF("shift=%d\n", shift);
713: fpu_shr(fp, shift - FP_NG);
714: DUMPFP("s2:shr ", fp);
715:
716: rup = fpu_round(fe, fp);
717: uint64_t m = (((uint64_t)fp->fp_mant[1]) << 32)
718: | (uint64_t)fp->fp_mant[2];
719: if (rup && m == (1ULL << effbits)) {
720: fp->fp_exp++;
721: shift--;
722: }
723: DUMPFP("s2:round", fp);
724:
725: fpu_shl(fp, shift);
726: DUMPFP("s2:shl ", fp);
727:
728: if (fp->fp_exp <= -EXP_BIAS - FRACBITS) {
729: PRINTF("s2:zero\n");
730: fp->fp_class = FPC_ZERO;
731: }
732: return;
733: }
734:
735: // 指定精度の正規化数で表現できる場合
736: DUMPFP("s3:start", fp);
737:
738: // [0] [1] [2]
739: // 8765432109876543210 0 0 0 0
740: // IMMMMMMMMMMMMMMMMMM m..m MMMMMMMMMMMMMMMMMMMMMMMMMMMMMMGR : origin
741: // IMMMMMMMMMMMMMMMMMMmmmmmmm : shifted
742: // 21098765432109876543210GR
743: shift = FP_NMANT - FRACBITS - 1/*Integer bit*/;
744: fpu_shr(fp, shift - FP_NG);
745: DUMPFP("s3:shr ", fp);
746:
747: rup = fpu_round(fe, fp);
748: uint64_t m = (((uint64_t)fp->fp_mant[1]) << 32)
749: | (uint64_t)fp->fp_mant[2];
750: if (rup && (m == (2ULL << FRACBITS))) {
751: fp->fp_exp++;
752: shift--;
753: }
754: DUMPFP("s3:round", fp);
755:
756: fpu_shl(fp, shift);
757: DUMPFP("s3:shl ", fp);
758: }
759:
760: void
761: fpu_round_chkinf(struct fpemu *fe, struct fpn *fp, int EXP_BIAS)
762: {
763: if (fp->fp_exp > EXP_BIAS) {
764: fe->fe_fpsr |= FPSR_OVFL | FPSR_INEX2;
765:
766: // 指定精度で表現できない大きい値の場合は
767: // 正で RN/RP か、負で RN/RM なら Inf、
768: // そうでなければ指定精度で表現できる最大値。
769: if (fp->fp_sign == 0) {
770: switch ((fe->fe_fpcr & FPCR_ROUND)) {
771: case FPCR_NEAR:
772: case FPCR_PINF:
773: fp->fp_class = FPC_INF;
774: break;
775: }
776: } else {
777: switch ((fe->fe_fpcr & FPCR_ROUND)) {
778: case FPCR_NEAR:
779: case FPCR_MINF:
780: fp->fp_class = FPC_INF;
781: break;
782: }
783: }
784: if (fp->fp_class != FPC_INF) {
785: // 最大値
786: fp->fp_exp = EXP_BIAS;
787: fp->fp_mant[0] = FP_2 -1;
788: switch ((fe->fe_fpcr & FPCR_PREC)) {
789: case FPCR_SNGL:
790: fp->fp_mant[1] = 0xf8000000;
791: fp->fp_mant[2] = 0;
792: break;
793: case FPCR_DBL:
794: fp->fp_mant[1] = 0xffffffff;
795: fp->fp_mant[2] = 0xc0000000;
796: break;
797: case FPCR_EXTD:
798: fp->fp_mant[1] = 0xffffffff;
799: fp->fp_mant[2] = 0xfff80000;
800: break;
801: }
802: }
803: }
804: }
805: #endif /* XM6i_FPE */
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.