5dea5ad2bc2cfc18a470f0db78dc21d898276398
1 /* vax_fpa.c - VAX f_, d_, g_floating instructions
3 Copyright (c) 1998-2008, Robert M Supnik
5 Permission is hereby granted, free of charge, to any person obtaining a
6 copy of this software and associated documentation files (the "Software"),
7 to deal in the Software without restriction, including without limitation
8 the rights to use, copy, modify, merge, publish, distribute, sublicense,
9 and/or sell copies of the Software, and to permit persons to whom the
10 Software is furnished to do so, subject to the following conditions:
12 The above copyright notice and this permission notice shall be included in
13 all copies or substantial portions of the Software.
15 THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
16 IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
17 FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
18 ROBERT M SUPNIK BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER
19 IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN
20 CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
22 Except as contained in this notice, the name of Robert M Supnik shall not be
23 used in advertising or otherwise to promote the sale, use or other dealings
24 in this Software without prior written authorization from Robert M Supnik.
26 28-May-08 RMS Inlined physical memory routines
27 16-May-06 RMS Fixed bug in 32b floating multiply routine
28 Fixed bug in 64b extended modulus routine
29 03-May-06 RMS Fixed POLYD, POLYG to clear R4, R5
30 Fixed POLYD, POLYG to set R3 correctly
31 Fixed POLYD, POLYG to not exit prematurely if arg = 0
32 Fixed POLYD, POLYG to do full 64b multiply
33 Fixed POLYF, POLYD, POLYG to remove truncation on add
34 Fixed POLYF, POLYD, POLYG to mask mul reslt to 31b/63b/63b
35 Fixed fp add routine to test for zero via fraction
36 to support "denormal" argument from POLYF, POLYD, POLYG
37 (all reported by Tim Stark)
38 27-Sep-05 RMS Fixed bug in 32b structure definitions (from Jason Stevens)
39 30-Sep-04 RMS Comment and formating changes based on vax_octa.c
40 18-Apr-04 RMS Moved format definitions to vax_defs.h
41 19-Jun-03 RMS Simplified add algorithm
42 16-May-03 RMS Fixed bug in floating to integer convert overflow
43 Fixed multiple bugs in EMODx
44 Integrated 32b only code
45 05-Jul-02 RMS Changed internal routine names for C library conflict
46 17-Apr-02 RMS Fixed bug in EDIV zero quotient
48 This module contains the instruction simulators for
50 - 64 bit arithmetic (ASHQ, EMUL, EDIV)
51 - single precision floating point
52 - double precision floating point, D and G format
61 extern jmp_buf save_env
;
63 #if defined (USE_INT64)
65 #define M64 0xFFFFFFFFFFFFFFFF /* 64b */
66 #define FD_FRACW (0xFFFF & ~(FD_EXP | FPSIGN))
67 #define FD_FRACL (FD_FRACW | 0xFFFF0000) /* f/d fraction */
68 #define G_FRACW (0xFFFF & ~(G_EXP | FPSIGN))
69 #define G_FRACL (G_FRACW | 0xFFFF0000) /* g fraction */
70 #define UNSCRAM(h,l) (((((t_uint64) (h)) << 48) & 0xFFFF000000000000) | \
71 ((((t_uint64) (h)) << 16) & 0x0000FFFF00000000) | \
72 ((((t_uint64) (l)) << 16) & 0x00000000FFFF0000) | \
73 ((((t_uint64) (l)) >> 16) & 0x000000000000FFFF))
74 #define CONCAT(h,l) ((((t_uint64) (h)) << 32) | ((uint32) (l)))
82 #define UF_NM 0x8000000000000000 /* normalized */
83 #define UF_FRND 0x0000008000000000 /* F round */
84 #define UF_DRND 0x0000000000000080 /* D round */
85 #define UF_GRND 0x0000000000000400 /* G round */
88 #define UF_V_FDLO (UF_V_FDHI - 32)
90 #define UF_V_GLO (UF_V_GHI - 32)
91 #define UF_GETFDHI(x) (int32) ((((x) >> (16 + UF_V_FDHI)) & FD_FRACW) | \
92 (((x) >> (UF_V_FDHI - 16)) & ~0xFFFF))
93 #define UF_GETFDLO(x) (int32) ((((x) >> (16 + UF_V_FDLO)) & 0xFFFF) | \
94 (((x) << (16 - UF_V_FDLO)) & ~0xFFFF))
95 #define UF_GETGHI(x) (int32) ((((x) >> (16 + UF_V_GHI)) & G_FRACW) | \
96 (((x) >> (UF_V_GHI - 16)) & ~0xFFFF))
97 #define UF_GETGLO(x) (int32) ((((x) >> (16 + UF_V_GLO)) & 0xFFFF) | \
98 (((x) << (16 - UF_V_GLO)) & ~0xFFFF))
100 void unpackf (int32 hi
, UFP
*a
);
101 void unpackd (int32 hi
, int32 lo
, UFP
*a
);
102 void unpackg (int32 hi
, int32 lo
, UFP
*a
);
104 int32
rpackfd (UFP
*a
, int32
*rh
);
105 int32
rpackg (UFP
*a
, int32
*rh
);
106 void vax_fadd (UFP
*a
, UFP
*b
);
107 void vax_fmul (UFP
*a
, UFP
*b
, t_bool qd
, int32 bias
, uint32 mhi
, uint32 mlo
);
108 void vax_fdiv (UFP
*b
, UFP
*a
, int32 prec
, int32 bias
);
109 void vax_fmod (UFP
*a
, int32 bias
, int32
*intgr
, int32
*flg
);
111 /* Quadword arithmetic shift
113 opnd[0] = shift count (cnt.rb)
114 opnd[1:2] = source (src.rq)
115 opnd[3:4] = destination (dst.wq)
118 int32
op_ashq (int32
*opnd
, int32
*rh
, int32
*flg
)
123 src
= CONCAT (opnd
[2], opnd
[1]); /* build src */
124 if (sc
& BSIGN
) { /* right shift? */
125 *flg
= 0; /* no ovflo */
126 sc
= 0x100 - sc
; /* |shift| */
127 if (sc
> 63) r
= (opnd
[2] & LSIGN
)? -1: 0; /* sc > 63? */
131 if (sc
> 63) { /* left shift */
132 r
= 0; /* sc > 63? */
133 *flg
= (src
!= 0); /* ovflo test */
136 r
= src
<< sc
; /* do shift */
137 *flg
= (src
!= (r
>> sc
)); /* ovflo test */
140 *rh
= (int32
) ((r
>> 32) & LMASK
); /* hi result */
141 return ((int32
) (r
& LMASK
)); /* lo result */
144 /* Extended multiply subroutine */
146 int32
op_emul (int32 mpy
, int32 mpc
, int32
*rh
)
152 *rh
= (int32
) ((lmpy
>> 32) & LMASK
);
153 return ((int32
) (lmpy
& LMASK
));
158 opnd[0] = divisor (non-zero)
162 int32
op_ediv (int32
*opnd
, int32
*rh
, int32
*flg
)
167 *flg
= CC_V
; /* assume error */
169 ldvr
= ((opnd
[0] & LSIGN
)? -opnd
[0]: opnd
[0]) & LMASK
; /* |divisor| */
170 ldvd
= CONCAT (opnd
[2], opnd
[1]); /* 64b dividend */
171 if (opnd
[2] & LSIGN
) ldvd
= -ldvd
; /* |dividend| */
172 if (((ldvd
>> 32) & LMASK
) >= ldvr
) return opnd
[1]; /* divide work? */
173 quo
= (int32
) (ldvd
/ ldvr
); /* do divide */
174 rem
= (int32
) (ldvd
% ldvr
);
175 if ((opnd
[0] ^ opnd
[2]) & LSIGN
) { /* result -? */
176 quo
= -quo
; /* negate */
177 if (quo
&& ((quo
& LSIGN
) == 0)) return opnd
[1]; /* right sign? */
179 else if (quo
& LSIGN
) return opnd
[1];
180 if (opnd
[2] & LSIGN
) rem
= -rem
; /* sign of rem */
181 *flg
= 0; /* no overflow */
182 *rh
= rem
& LMASK
; /* set rem */
183 return (quo
& LMASK
); /* return quo */
186 /* Compare floating */
188 int32
op_cmpfd (int32 h1
, int32 l1
, int32 h2
, int32 l2
)
192 if ((h1
& FD_EXP
) == 0) {
193 if (h1
& FPSIGN
) RSVD_OPND_FAULT
;
196 if ((h2
& FD_EXP
) == 0) {
197 if (h2
& FPSIGN
) RSVD_OPND_FAULT
;
200 if ((h1
^ h2
) & FPSIGN
) return ((h1
& FPSIGN
)? CC_N
: 0);
201 n1
= UNSCRAM (h1
, l1
);
202 n2
= UNSCRAM (h2
, l2
);
203 if (n1
== n2
) return CC_Z
;
204 return (((n1
< n2
) ^ ((h1
& FPSIGN
) != 0))? CC_N
: 0);
207 int32
op_cmpg (int32 h1
, int32 l1
, int32 h2
, int32 l2
)
211 if ((h1
& G_EXP
) == 0) {
212 if (h1
& FPSIGN
) RSVD_OPND_FAULT
;
215 if ((h2
& G_EXP
) == 0) {
216 if (h2
& FPSIGN
) RSVD_OPND_FAULT
;
219 if ((h1
^ h2
) & FPSIGN
) return ((h1
& FPSIGN
)? CC_N
: 0);
220 n1
= UNSCRAM (h1
, l1
);
221 n2
= UNSCRAM (h2
, l2
);
222 if (n1
== n2
) return CC_Z
;
223 return (((n1
< n2
) ^ ((h1
& FPSIGN
) != 0))? CC_N
: 0);
226 /* Integer to floating convert */
228 int32
op_cvtifdg (int32 val
, int32
*rh
, int32 opc
)
241 a
.exp
= 32 + ((opc
& 0x100)? G_BIAS
: FD_BIAS
);
242 a
.frac
= ((t_uint64
) val
) << (UF_V_NM
- 31);
244 if (opc
& 0x100) return rpackg (&a
, rh
);
245 return rpackfd (&a
, rh
);
248 /* Floating to integer convert */
250 int32
op_cvtfdgi (int32
*opnd
, int32
*flg
, int32 opc
)
253 int32 lnt
= opc
& 03;
255 static t_uint64 maxv
[4] = { 0x7F, 0x7FFF, 0x7FFFFFFF, 0x7FFFFFFF };
259 unpackg (opnd
[0], opnd
[1], &a
);
260 ubexp
= a
.exp
- G_BIAS
;
263 if (opc
& 0x20) unpackd (opnd
[0], opnd
[1], &a
);
264 else unpackf (opnd
[0], &a
);
265 ubexp
= a
.exp
- FD_BIAS
;
267 if ((a
.exp
== 0) || (ubexp
< 0)) return 0;
268 if (ubexp
<= UF_V_NM
) {
269 a
.frac
= a
.frac
>> (UF_V_NM
- ubexp
); /* leave rnd bit */
270 if ((opc
& 03) == 03) a
.frac
= a
.frac
+ 1; /* if CVTR, round */
271 a
.frac
= a
.frac
>> 1; /* now justified */
272 if (a
.frac
> (maxv
[lnt
] + (a
.sign
? 1: 0))) *flg
= CC_V
;
275 *flg
= CC_V
; /* set overflow */
276 if (ubexp
> (UF_V_NM
+ 32)) return 0;
277 a
.frac
= a
.frac
<< (ubexp
- UF_V_NM
- 1); /* no rnd bit */
279 return ((int32
) ((a
.sign
? (a
.frac
^ LMASK
) + 1: a
.frac
) & LMASK
));
282 /* Extended modularize
284 One of three floating point instructions dropped from the architecture,
285 EMOD presents two sets of complications. First, it requires an extended
286 fraction multiply, with precise (and unusual) truncation conditions.
287 Second, it has two write operands, a dubious distinction it shares
291 int32
op_emodf (int32
*opnd
, int32
*intgr
, int32
*flg
)
295 unpackf (opnd
[0], &a
); /* unpack operands */
296 unpackf (opnd
[2], &b
);
297 a
.frac
= a
.frac
| (((t_uint64
) opnd
[1]) << 32); /* extend src1 */
298 vax_fmul (&a
, &b
, 0, FD_BIAS
, 0, LMASK
); /* multiply */
299 vax_fmod (&a
, FD_BIAS
, intgr
, flg
); /* sep int & frac */
300 return rpackfd (&a
, NULL
); /* return frac */
303 int32
op_emodd (int32
*opnd
, int32
*flo
, int32
*intgr
, int32
*flg
)
307 unpackd (opnd
[0], opnd
[1], &a
); /* unpack operands */
308 unpackd (opnd
[3], opnd
[4], &b
);
309 a
.frac
= a
.frac
| opnd
[2]; /* extend src1 */
310 vax_fmul (&a
, &b
, 1, FD_BIAS
, 0, 0); /* multiply */
311 vax_fmod (&a
, FD_BIAS
, intgr
, flg
); /* sep int & frac */
312 return rpackfd (&a
, flo
); /* return frac */
315 int32
op_emodg (int32
*opnd
, int32
*flo
, int32
*intgr
, int32
*flg
)
319 unpackg (opnd
[0], opnd
[1], &a
); /* unpack operands */
320 unpackg (opnd
[3], opnd
[4], &b
);
321 a
.frac
= a
.frac
| (opnd
[2] >> 5); /* extend src1 */
322 vax_fmul (&a
, &b
, 1, G_BIAS
, 0, 0); /* multiply */
323 vax_fmod (&a
, G_BIAS
, intgr
, flg
); /* sep int & frac */
324 return rpackg (&a
, flo
); /* return frac */
327 /* Unpacked floating point routines */
329 void vax_fadd (UFP
*a
, UFP
*b
)
334 if (a
->frac
== 0) { /* s1 = 0? */
338 if (b
->frac
== 0) return; /* s2 = 0? */
339 if ((a
->exp
< b
->exp
) || /* |s1| < |s2|? swap */
340 ((a
->exp
== b
->exp
) && (a
->frac
< b
->frac
))) {
345 ediff
= a
->exp
- b
->exp
; /* exp diff */
346 if (a
->sign
^ b
->sign
) { /* eff sub? */
347 if (ediff
) { /* exp diff? */
348 if (ediff
> 63) b
->frac
= M64
; /* retain sticky */
349 else b
->frac
= ((-((t_int64
) b
->frac
) >> ediff
) | /* denormalize */
350 (M64
<< (64 - ediff
))); /* preserve sign */
351 a
->frac
= a
->frac
+ b
->frac
; /* add frac */
353 else a
->frac
= a
->frac
- b
->frac
; /* sub frac */
354 norm (a
); /* normalize */
357 if (ediff
> 63) b
->frac
= 0; /* add */
358 else if (ediff
) b
->frac
= b
->frac
>> ediff
; /* denormalize */
359 a
->frac
= a
->frac
+ b
->frac
; /* add frac */
360 if (a
->frac
< b
->frac
) { /* chk for carry */
361 a
->frac
= UF_NM
| (a
->frac
>> 1); /* shift in carry */
362 a
->exp
= a
->exp
+ 1; /* skip norm */
368 /* Floating multiply - 64b * 64b with cross products */
370 void vax_fmul (UFP
*a
, UFP
*b
, t_bool qd
, int32 bias
, uint32 mhi
, uint32 mlo
)
372 t_uint64 ah
, bh
, al
, bl
, rhi
, rlo
, rmid1
, rmid2
;
373 t_uint64 mask
= (((t_uint64
) mhi
) << 32) | ((t_uint64
) mlo
);
375 if ((a
->exp
== 0) || (b
->exp
== 0)) { /* zero argument? */
376 a
->frac
= a
->sign
= a
->exp
= 0; /* result is zero */
379 a
->sign
= a
->sign
^ b
->sign
; /* sign of result */
380 a
->exp
= a
->exp
+ b
->exp
- bias
; /* add exponents */
381 ah
= (a
->frac
>> 32) & LMASK
; /* split operands */
382 bh
= (b
->frac
>> 32) & LMASK
; /* into 32b chunks */
383 rhi
= ah
* bh
; /* high result */
384 if (qd
) { /* 64b needed? */
385 al
= a
->frac
& LMASK
;
386 bl
= b
->frac
& LMASK
;
390 rhi
= rhi
+ ((rmid1
>> 32) & LMASK
) + ((rmid2
>> 32) & LMASK
);
391 rmid1
= rlo
+ (rmid1
<< 32); /* add mid1 to lo */
392 if (rmid1
< rlo
) rhi
= rhi
+ 1; /* carry? incr hi */
393 rmid2
= rmid1
+ (rmid2
<< 32); /* add mid2 to lo */
394 if (rmid2
< rmid1
) rhi
= rhi
+ 1; /* carry? incr hi */
396 a
->frac
= rhi
& ~mask
;
397 norm (a
); /* normalize */
401 /* Floating modulus - there are three cases
403 exp <= bias - integer is 0, fraction is input,
405 bias < exp <= bias+64 - separate integer and fraction,
406 integer overflow may occur
407 bias+64 < exp - result is integer, fraction is 0
411 void vax_fmod (UFP
*a
, int32 bias
, int32
*intgr
, int32
*flg
)
413 if (a
->exp
<= bias
) *intgr
= *flg
= 0; /* 0 or <1? int = 0 */
414 else if (a
->exp
<= (bias
+ 64)) { /* in range [1,64]? */
415 *intgr
= (int32
) (a
->frac
>> (64 - (a
->exp
- bias
)));
416 if ((a
->exp
> (bias
+ 32)) || /* test ovflo */
417 ((a
->exp
== (bias
+ 32)) &&
418 (((uint32
) *intgr
) > (a
->sign
? 0x80000000: 0x7FFFFFFF))))
421 if (a
->sign
) *intgr
= -*intgr
; /* -? comp int */
422 if (a
->exp
== (bias
+ 64)) a
->frac
= 0; /* special case 64 */
423 else a
->frac
= a
->frac
<< (a
->exp
- bias
);
427 *intgr
= 0; /* out of range */
428 a
->frac
= a
->sign
= a
->exp
= 0; /* result 0 */
429 *flg
= CC_V
; /* overflow */
431 norm (a
); /* normalize */
436 Needs to develop at least one rounding bit. Since the first
437 divide step can fail, caller should specify 2 more bits than
438 the precision of the fraction.
441 void vax_fdiv (UFP
*a
, UFP
*b
, int32 prec
, int32 bias
)
446 if (a
->exp
== 0) FLT_DZRO_FAULT
; /* divr = 0? */
447 if (b
->exp
== 0) return; /* divd = 0? */
448 b
->sign
= b
->sign
^ a
->sign
; /* result sign */
449 b
->exp
= b
->exp
- a
->exp
+ bias
+ 1; /* unbiased exp */
450 a
->frac
= a
->frac
>> 1; /* allow 1 bit left */
451 b
->frac
= b
->frac
>> 1;
452 for (i
= 0; (i
< prec
) && b
->frac
; i
++) { /* divide loop */
453 quo
= quo
<< 1; /* shift quo */
454 if (b
->frac
>= a
->frac
) { /* div step ok? */
455 b
->frac
= b
->frac
- a
->frac
; /* subtract */
456 quo
= quo
+ 1; /* quo bit = 1 */
458 b
->frac
= b
->frac
<< 1; /* shift divd */
460 b
->frac
= quo
<< (UF_V_NM
- i
+ 1); /* shift quo */
461 norm (b
); /* normalize */
465 /* Support routines */
467 void unpackf (int32 hi
, UFP
*r
)
469 r
->sign
= hi
& FPSIGN
; /* get sign */
470 r
->exp
= FD_GETEXP (hi
); /* get exponent */
471 if (r
->exp
== 0) { /* exp = 0? */
472 if (r
->sign
) RSVD_OPND_FAULT
; /* if -, rsvd op */
473 r
->frac
= 0; /* else 0 */
476 hi
= (((hi
& FD_FRACW
) | FD_HB
) << 16) | ((hi
>> 16) & 0xFFFF);
477 r
->frac
= ((t_uint64
) hi
) << (32 + UF_V_FDLO
);
481 void unpackd (int32 hi
, int32 lo
, UFP
*r
)
483 r
->sign
= hi
& FPSIGN
; /* get sign */
484 r
->exp
= FD_GETEXP (hi
); /* get exponent */
485 if (r
->exp
== 0) { /* exp = 0? */
486 if (r
->sign
) RSVD_OPND_FAULT
; /* if -, rsvd op */
487 r
->frac
= 0; /* else 0 */
490 hi
= (hi
& FD_FRACL
) | FD_HB
; /* canonical form */
491 r
->frac
= UNSCRAM (hi
, lo
) << UF_V_FDLO
; /* guard bits */
495 void unpackg (int32 hi
, int32 lo
, UFP
*r
)
497 r
->sign
= hi
& FPSIGN
; /* get sign */
498 r
->exp
= G_GETEXP (hi
); /* get exponent */
499 if (r
->exp
== 0) { /* exp = 0? */
500 if (r
->sign
) RSVD_OPND_FAULT
; /* if -, rsvd op */
501 r
->frac
= 0; /* else 0 */
504 hi
= (hi
& G_FRACL
) | G_HB
; /* canonical form */
505 r
->frac
= UNSCRAM (hi
, lo
) << UF_V_GLO
; /* guard bits */
512 static t_uint64 normmask
[5] = {
513 0xc000000000000000, 0xf000000000000000, 0xff00000000000000,
514 0xffff000000000000, 0xffffffff00000000
516 static int32 normtab
[6] = { 1, 2, 4, 8, 16, 32};
518 if (r
->frac
== 0) { /* if fraction = 0 */
519 r
->sign
= r
->exp
= 0; /* result is 0 */
522 while ((r
->frac
& UF_NM
) == 0) { /* normalized? */
523 for (i
= 0; i
< 5; i
++) { /* find first 1 */
524 if (r
->frac
& normmask
[i
]) break;
526 r
->frac
= r
->frac
<< normtab
[i
]; /* shift frac */
527 r
->exp
= r
->exp
- normtab
[i
]; /* decr exp */
532 int32
rpackfd (UFP
*r
, int32
*rh
)
534 if (rh
) *rh
= 0; /* assume 0 */
535 if (r
->frac
== 0) return 0; /* result 0? */
536 r
->frac
= r
->frac
+ (rh
? UF_DRND
: UF_FRND
); /* round */
537 if ((r
->frac
& UF_NM
) == 0) { /* carry out? */
538 r
->frac
= r
->frac
>> 1; /* renormalize */
541 if (r
->exp
> (int32
) FD_M_EXP
) FLT_OVFL_FAULT
; /* ovflo? fault */
542 if (r
->exp
<= 0) { /* underflow? */
543 if (PSL
& PSW_FU
) FLT_UNFL_FAULT
; /* fault if fu */
544 return 0; /* else 0 */
546 if (rh
) *rh
= UF_GETFDLO (r
->frac
); /* get low */
547 return r
->sign
| (r
->exp
<< FD_V_EXP
) | UF_GETFDHI (r
->frac
);
550 int32
rpackg (UFP
*r
, int32
*rh
)
552 *rh
= 0; /* assume 0 */
553 if (r
->frac
== 0) return 0; /* result 0? */
554 r
->frac
= r
->frac
+ UF_GRND
; /* round */
555 if ((r
->frac
& UF_NM
) == 0) { /* carry out? */
556 r
->frac
= r
->frac
>> 1; /* renormalize */
559 if (r
->exp
> (int32
) G_M_EXP
) FLT_OVFL_FAULT
; /* ovflo? fault */
560 if (r
->exp
<= 0) { /* underflow? */
561 if (PSL
& PSW_FU
) FLT_UNFL_FAULT
; /* fault if fu */
562 return 0; /* else 0 */
564 if (rh
) *rh
= UF_GETGLO (r
->frac
); /* get low */
565 return r
->sign
| (r
->exp
<< G_V_EXP
) | UF_GETGHI (r
->frac
);
570 #define WORDSWAP(x) ((((x) & WMASK) << 16) | (((x) >> 16) & WMASK))
583 #define UF_NM_H 0x80000000 /* normalized */
584 #define UF_FRND_H 0x00000080 /* F round */
585 #define UF_FRND_L 0x00000000
586 #define UF_DRND_H 0x00000000 /* D round */
587 #define UF_DRND_L 0x00000080
588 #define UF_GRND_H 0x00000000 /* G round */
589 #define UF_GRND_L 0x00000400
592 void unpackf (uint32 hi
, UFP
*a
);
593 void unpackd (uint32 hi
, uint32 lo
, UFP
*a
);
594 void unpackg (uint32 hi
, uint32 lo
, UFP
*a
);
596 int32
rpackfd (UFP
*a
, int32
*rh
);
597 int32
rpackg (UFP
*a
, int32
*rh
);
598 void vax_fadd (UFP
*a
, UFP
*b
);
599 void vax_fmul (UFP
*a
, UFP
*b
, t_bool qd
, int32 bias
, uint32 mhi
, uint32 mlo
);
600 void vax_fmod (UFP
*a
, int32 bias
, int32
*intgr
, int32
*flg
);
601 void vax_fdiv (UFP
*b
, UFP
*a
, int32 prec
, int32 bias
);
602 void dp_add (UDP
*a
, UDP
*b
);
603 void dp_inc (UDP
*a
);
604 void dp_sub (UDP
*a
, UDP
*b
);
605 void dp_imul (uint32 a
, uint32 b
, UDP
*r
);
606 void dp_lsh (UDP
*a
, uint32 sc
);
607 void dp_rsh (UDP
*a
, uint32 sc
);
608 void dp_rsh_s (UDP
*a
, uint32 sc
, uint32 neg
);
609 void dp_neg (UDP
*a
);
610 int32
dp_cmp (UDP
*a
, UDP
*b
);
612 /* Quadword arithmetic shift
614 opnd[0] = shift count (cnt.rb)
615 opnd[1:2] = source (src.rq)
616 opnd[3:4] = destination (dst.wq)
619 int32
op_ashq (int32
*opnd
, int32
*rh
, int32
*flg
)
624 r
.lo
= opnd
[1]; /* get source */
626 *flg
= 0; /* assume no ovflo */
627 if (sc
& BSIGN
) /* right shift? */
628 dp_rsh_s (&r
, 0x100 - sc
, r
.hi
& LSIGN
); /* signed right */
630 dp_lsh (&r
, sc
); /* left shift */
631 sr
= r
; /* copy result */
632 dp_rsh_s (&sr
, sc
, sr
.hi
& LSIGN
); /* signed right */
633 if ((sr
.hi
!= ((uint32
) opnd
[2])) || /* reshift != orig? */
634 (sr
.lo
!= ((uint32
) opnd
[1]))) *flg
= 1; /* overflow */
636 *rh
= r
.hi
; /* hi result */
637 return r
.lo
; /* lo result */
640 /* Extended multiply subroutine */
642 int32
op_emul (int32 mpy
, int32 mpc
, int32
*rh
)
645 int32 sign
= mpy
^ mpc
; /* sign of result */
647 if (mpy
& LSIGN
) mpy
= -mpy
; /* abs value */
648 if (mpc
& LSIGN
) mpc
= -mpc
;
649 dp_imul (mpy
& LMASK
, mpc
& LMASK
, &r
); /* 32b * 32b -> 64b */
650 if (sign
& LSIGN
) dp_neg (&r
); /* negative result? */
657 opnd[0] = divisor (non-zero)
661 int32
op_ediv (int32
*opnd
, int32
*rh
, int32
*flg
)
666 dvr
= opnd
[0]; /* get divisor */
667 dvd
.lo
= opnd
[1]; /* get dividend */
669 *flg
= CC_V
; /* assume error */
671 if (dvd
.hi
& LSIGN
) dp_neg (&dvd
); /* |dividend| */
672 if (dvr
& LSIGN
) dvr
= NEG (dvr
); /* |divisor| */
673 if (dvd
.hi
>= dvr
) return opnd
[1]; /* divide work? */
674 for (i
= quo
= 0; i
< 32; i
++) { /* 32 iterations */
675 quo
= quo
<< 1; /* shift quotient */
676 dp_lsh (&dvd
, 1); /* shift dividend */
677 if (dvd
.hi
>= dvr
) { /* step work? */
678 dvd
.hi
= (dvd
.hi
- dvr
) & LMASK
; /* subtract dvr */
682 if ((opnd
[0] ^ opnd
[2]) & LSIGN
) { /* result -? */
683 quo
= NEG (quo
); /* negate */
684 if (quo
&& ((quo
& LSIGN
) == 0)) return opnd
[1]; /* right sign? */
686 else if (quo
& LSIGN
) return opnd
[1];
687 if (opnd
[2] & LSIGN
) *rh
= NEG (dvd
.hi
); /* sign of rem */
689 *flg
= 0; /* no overflow */
690 return quo
; /* return quo */
693 /* Compare floating */
695 int32
op_cmpfd (int32 h1
, int32 l1
, int32 h2
, int32 l2
)
700 unpackd (h1
, l1
, &a
);
701 unpackd (h2
, l2
, &b
);
702 if (a
.sign
!= b
.sign
) return (a
.sign
? CC_N
: 0);
704 if (r
== 0) r
= dp_cmp (&a
.frac
, &b
.frac
);
705 if (r
< 0) return (a
.sign
? 0: CC_N
);
706 if (r
> 0) return (a
.sign
? CC_N
: 0);
710 int32
op_cmpg (int32 h1
, int32 l1
, int32 h2
, int32 l2
)
715 unpackg (h1
, l1
, &a
);
716 unpackg (h2
, l2
, &b
);
717 if (a
.sign
!= b
.sign
) return (a
.sign
? CC_N
: 0);
719 if (r
== 0) r
= dp_cmp (&a
.frac
, &b
.frac
);
720 if (r
< 0) return (a
.sign
? 0: CC_N
);
721 if (r
> 0) return (a
.sign
? CC_N
: 0);
725 /* Integer to floating convert */
727 int32
op_cvtifdg (int32 val
, int32
*rh
, int32 opc
)
731 if (val
== 0) { /* zero? */
732 if (rh
) *rh
= 0; /* return true 0 */
735 if (val
< 0) { /* negative? */
736 a
.sign
= FPSIGN
; /* sign = - */
739 else a
.sign
= 0; /* else sign = + */
740 a
.exp
= 32 + ((opc
& 0x100)? G_BIAS
: FD_BIAS
); /* initial exp */
741 a
.frac
.hi
= val
& LMASK
; /* fraction */
743 norm (&a
); /* normalize */
744 if (opc
& 0x100) return rpackg (&a
, rh
); /* pack and return */
745 return rpackfd (&a
, rh
);
748 /* Floating to integer convert */
750 int32
op_cvtfdgi (int32
*opnd
, int32
*flg
, int32 opc
)
753 int32 lnt
= opc
& 03;
755 static uint32 maxv
[4] = { 0x7F, 0x7FFF, 0x7FFFFFFF, 0x7FFFFFFF };
758 if (opc
& 0x100) { /* G? */
759 unpackg (opnd
[0], opnd
[1], &a
); /* unpack */
760 ubexp
= a
.exp
- G_BIAS
; /* unbiased exp */
763 if (opc
& 0x20) unpackd (opnd
[0], opnd
[1], &a
); /* F or D */
764 else unpackf (opnd
[0], &a
); /* unpack */
765 ubexp
= a
.exp
- FD_BIAS
; /* unbiased exp */
767 if ((a
.exp
== 0) || (ubexp
< 0)) return 0; /* true zero or frac? */
768 if (ubexp
<= UF_V_NM
) { /* exp in range? */
769 dp_rsh (&a
.frac
, UF_V_NM
- ubexp
); /* leave rnd bit */
770 if (lnt
== 03) dp_inc (&a
.frac
); /* if CVTR, round */
771 dp_rsh (&a
.frac
, 1); /* now justified */
772 if ((a
.frac
.hi
!= 0) ||
773 (a
.frac
.lo
> (maxv
[lnt
] + (a
.sign
? 1: 0)))) *flg
= CC_V
;
776 *flg
= CC_V
; /* always ovflo */
777 if (ubexp
> (UF_V_NM
+ 32)) return 0; /* in ext range? */
778 dp_lsh (&a
.frac
, ubexp
- UF_V_NM
- 1); /* no rnd bit */
780 return (a
.sign
? NEG (a
.frac
.lo
): a
.frac
.lo
); /* return lo frac */
783 /* Extended modularize
785 One of three floating point instructions dropped from the architecture,
786 EMOD presents two sets of complications. First, it requires an extended
787 fraction multiply, with precise (and unusual) truncation conditions.
788 Second, it has two write operands, a dubious distinction it shares
792 int32
op_emodf (int32
*opnd
, int32
*intgr
, int32
*flg
)
796 unpackf (opnd
[0], &a
); /* unpack operands */
797 unpackf (opnd
[2], &b
);
798 a
.frac
.hi
= a
.frac
.hi
| opnd
[1]; /* extend src1 */
799 vax_fmul (&a
, &b
, 0, FD_BIAS
, 0, LMASK
); /* multiply */
800 vax_fmod (&a
, FD_BIAS
, intgr
, flg
); /* sep int & frac */
801 return rpackfd (&a
, NULL
); /* return frac */
804 int32
op_emodd (int32
*opnd
, int32
*flo
, int32
*intgr
, int32
*flg
)
808 unpackd (opnd
[0], opnd
[1], &a
); /* unpack operands */
809 unpackd (opnd
[3], opnd
[4], &b
);
810 a
.frac
.lo
= a
.frac
.lo
| opnd
[2]; /* extend src1 */
811 vax_fmul (&a
, &b
, 1, FD_BIAS
, 0, 0); /* multiply */
812 vax_fmod (&a
, FD_BIAS
, intgr
, flg
); /* sep int & frac */
813 return rpackfd (&a
, flo
); /* return frac */
816 int32
op_emodg (int32
*opnd
, int32
*flo
, int32
*intgr
, int32
*flg
)
820 unpackg (opnd
[0], opnd
[1], &a
); /* unpack operands */
821 unpackg (opnd
[3], opnd
[4], &b
);
822 a
.frac
.lo
= a
.frac
.lo
| (opnd
[2] >> 5); /* extend src1 */
823 vax_fmul (&a
, &b
, 1, G_BIAS
, 0, 0); /* multiply */
824 vax_fmod (&a
, G_BIAS
, intgr
, flg
); /* sep int & frac */
825 return rpackg (&a
, flo
); /* return frac */
828 /* Unpacked floating point routines */
832 void vax_fadd (UFP
*a
, UFP
*b
)
837 if ((a
->frac
.hi
== 0) && (a
->frac
.lo
== 0)) { /* s1 = 0? */
841 if ((b
->frac
.hi
== 0) && (b
->frac
.lo
== 0)) return; /* s2 = 0? */
842 if ((a
->exp
< b
->exp
) || /* |s1| < |s2|? swap */
843 ((a
->exp
== b
->exp
) && (dp_cmp (&a
->frac
, &b
->frac
) < 0))) {
848 ediff
= a
->exp
- b
->exp
; /* exp diff */
849 if (a
->sign
^ b
->sign
) { /* eff sub? */
850 if (ediff
) { /* exp diff? */
851 dp_neg (&b
->frac
); /* negate fraction */
852 dp_rsh_s (&b
->frac
, ediff
, 1); /* signed right */
853 dp_add (&a
->frac
, &b
->frac
); /* "add" frac */
855 else dp_sub (&a
->frac
, &b
->frac
); /* a >= b */
856 norm (a
); /* normalize */
859 if (ediff
) dp_rsh (&b
->frac
, ediff
); /* add, denormalize */
860 dp_add (&a
->frac
, &b
->frac
); /* add frac */
861 if (dp_cmp (&a
->frac
, &b
->frac
) < 0) { /* chk for carry */
862 dp_rsh (&a
->frac
, 1); /* renormalize */
863 a
->frac
.hi
= a
->frac
.hi
| UF_NM_H
; /* add norm bit */
864 a
->exp
= a
->exp
+ 1; /* skip norm */
870 /* Floating multiply - 64b * 64b with cross products */
872 void vax_fmul (UFP
*a
, UFP
*b
, t_bool qd
, int32 bias
, uint32 mhi
, uint32 mlo
)
874 UDP rhi
, rlo
, rmid1
, rmid2
;
876 if ((a
->exp
== 0) || (b
->exp
== 0)) { /* zero argument? */
877 a
->frac
.hi
= a
->frac
.lo
= 0; /* result is zero */
878 a
->sign
= a
->exp
= 0;
881 a
->sign
= a
->sign
^ b
->sign
; /* sign of result */
882 a
->exp
= a
->exp
+ b
->exp
- bias
; /* add exponents */
883 dp_imul (a
->frac
.hi
, b
->frac
.hi
, &rhi
); /* high result */
884 if (qd
) { /* 64b needed? */
885 dp_imul (a
->frac
.hi
, b
->frac
.lo
, &rmid1
); /* cross products */
886 dp_imul (a
->frac
.lo
, b
->frac
.hi
, &rmid2
);
887 dp_imul (a
->frac
.lo
, b
->frac
.lo
, &rlo
); /* low result */
888 rhi
.lo
= (rhi
.lo
+ rmid1
.hi
) & LMASK
; /* add hi cross */
889 if (rhi
.lo
< rmid1
.hi
) /* to low high res */
890 rhi
.hi
= (rhi
.hi
+ 1) & LMASK
;
891 rhi
.lo
= (rhi
.lo
+ rmid2
.hi
) & LMASK
;
892 if (rhi
.lo
< rmid2
.hi
)
893 rhi
.hi
= (rhi
.hi
+ 1) & LMASK
;
894 rlo
.hi
= (rlo
.hi
+ rmid1
.lo
) & LMASK
; /* add mid1 to low res */
895 if (rlo
.hi
< rmid1
.lo
) dp_inc (&rhi
); /* carry? incr high res */
896 rlo
.hi
= (rlo
.hi
+ rmid2
.lo
) & LMASK
; /* add mid2 to low res */
897 if (rlo
.hi
< rmid2
.lo
) dp_inc (&rhi
); /* carry? incr high res */
899 a
->frac
.hi
= rhi
.hi
& ~mhi
; /* mask fraction */
900 a
->frac
.lo
= rhi
.lo
& ~mlo
;
901 norm (a
); /* normalize */
905 /* Floating modulus - there are three cases
907 exp <= bias - integer is 0, fraction is input,
909 bias < exp <= bias+64 - separate integer and fraction,
910 integer overflow may occur
911 bias+64 < exp - result is integer, fraction is 0
915 void vax_fmod (UFP
*a
, int32 bias
, int32
*intgr
, int32
*flg
)
919 if (a
->exp
<= bias
) *intgr
= *flg
= 0; /* 0 or <1? int = 0 */
920 else if (a
->exp
<= (bias
+ 64)) { /* in range [1,64]? */
922 dp_rsh (&ifr
, 64 - (a
->exp
- bias
)); /* separate integer */
923 if ((a
->exp
> (bias
+ 32)) || /* test ovflo */
924 ((a
->exp
== (bias
+ 32)) &&
925 (ifr
.lo
> (a
->sign
? 0x80000000: 0x7FFFFFFF))))
929 if (a
->sign
) *intgr
= -*intgr
; /* -? comp int */
930 dp_lsh (&a
->frac
, a
->exp
- bias
); /* excise integer */
934 *intgr
= 0; /* out of range */
935 a
->frac
.hi
= a
->frac
.lo
= a
->sign
= a
->exp
= 0; /* result 0 */
936 *flg
= CC_V
; /* overflow */
938 norm (a
); /* normalize */
943 Needs to develop at least one rounding bit. Since the first
944 divide step can fail, caller should specify 2 more bits than
945 the precision of the fraction.
948 void vax_fdiv (UFP
*a
, UFP
*b
, int32 prec
, int32 bias
)
953 if (a
->exp
== 0) FLT_DZRO_FAULT
; /* divr = 0? */
954 if (b
->exp
== 0) return; /* divd = 0? */
955 b
->sign
= b
->sign
^ a
->sign
; /* result sign */
956 b
->exp
= b
->exp
- a
->exp
+ bias
+ 1; /* unbiased exp */
957 dp_rsh (&a
->frac
, 1); /* allow 1 bit left */
958 dp_rsh (&b
->frac
, 1);
959 for (i
= 0; i
< prec
; i
++) { /* divide loop */
960 dp_lsh (&quo
, 1); /* shift quo */
961 if (dp_cmp (&b
->frac
, &a
->frac
) >= 0) { /* div step ok? */
962 dp_sub (&b
->frac
, &a
->frac
); /* subtract */
963 quo
.lo
= quo
.lo
+ 1; /* quo bit = 1 */
965 dp_lsh (&b
->frac
, 1); /* shift divd */
967 dp_lsh (&quo
, UF_V_NM
- prec
+ 1); /* put in position */
969 norm (b
); /* normalize */
973 /* Double precision integer routines */
975 int32
dp_cmp (UDP
*a
, UDP
*b
)
977 if (a
->hi
< b
->hi
) return -1; /* compare hi */
978 if (a
->hi
> b
->hi
) return +1;
979 if (a
->lo
< b
->lo
) return -1; /* hi =, compare lo */
980 if (a
->lo
> b
->lo
) return +1;
981 return 0; /* hi, lo equal */
984 void dp_add (UDP
*a
, UDP
*b
)
986 a
->lo
= (a
->lo
+ b
->lo
) & LMASK
; /* add lo */
987 if (a
->lo
< b
->lo
) a
->hi
= a
->hi
+ 1; /* carry? */
988 a
->hi
= (a
->hi
+ b
->hi
) & LMASK
; /* add hi */
994 a
->lo
= (a
->lo
+ 1) & LMASK
; /* inc lo */
995 if (a
->lo
== 0) a
->hi
= (a
->hi
+ 1) & LMASK
; /* carry? inc hi */
999 void dp_sub (UDP
*a
, UDP
*b
)
1001 if (a
->lo
< b
->lo
) a
->hi
= a
->hi
- 1; /* borrow? decr hi */
1002 a
->lo
= (a
->lo
- b
->lo
) & LMASK
; /* sub lo */
1003 a
->hi
= (a
->hi
- b
->hi
) & LMASK
; /* sub hi */
1007 void dp_lsh (UDP
*r
, uint32 sc
)
1009 if (sc
> 63) r
->hi
= r
->lo
= 0; /* > 63? result 0 */
1010 else if (sc
> 31) { /* [32,63]? */
1011 r
->hi
= (r
->lo
<< (sc
- 32)) & LMASK
;
1015 r
->hi
= ((r
->hi
<< sc
) | (r
->lo
>> (32 - sc
))) & LMASK
;
1016 r
->lo
= (r
->lo
<< sc
) & LMASK
;
1021 void dp_rsh (UDP
*r
, uint32 sc
)
1023 if (sc
> 63) r
->hi
= r
->lo
= 0; /* > 63? result 0 */
1024 else if (sc
> 31) { /* [32,63]? */
1025 r
->lo
= (r
->hi
>> (sc
- 32)) & LMASK
;
1029 r
->lo
= ((r
->lo
>> sc
) | (r
->hi
<< (32 - sc
))) & LMASK
;
1030 r
->hi
= (r
->hi
>> sc
) & LMASK
;
1035 void dp_rsh_s (UDP
*r
, uint32 sc
, uint32 neg
)
1037 dp_rsh (r
, sc
); /* do unsigned right */
1038 if (neg
&& sc
) { /* negative? */
1039 if (sc
> 63) r
->hi
= r
->lo
= LMASK
; /* > 63? result -1 */
1041 UDP ones
= { LMASK
, LMASK
};
1042 dp_lsh (&ones
, 64 - sc
); /* shift ones */
1043 r
->hi
= r
->hi
| ones
.hi
; /* or into result */
1044 r
->lo
= r
->lo
| ones
.lo
;
1050 void dp_imul (uint32 a
, uint32 b
, UDP
*r
)
1052 uint32 ah
, bh
, al
, bl
, rhi
, rlo
, rmid1
, rmid2
;
1054 if ((a
== 0) || (b
== 0)) { /* zero argument? */
1055 r
->hi
= r
->lo
= 0; /* result is zero */
1058 ah
= (a
>> 16) & WMASK
; /* split operands */
1059 bh
= (b
>> 16) & WMASK
; /* into 16b chunks */
1062 rhi
= ah
* bh
; /* high result */
1066 rhi
= rhi
+ ((rmid1
>> 16) & WMASK
) + ((rmid2
>> 16) & WMASK
);
1067 rmid1
= (rlo
+ (rmid1
<< 16)) & LMASK
; /* add mid1 to lo */
1068 if (rmid1
< rlo
) rhi
= rhi
+ 1; /* carry? incr hi */
1069 rmid2
= (rmid1
+ (rmid2
<< 16)) & LMASK
; /* add mid2 to to */
1070 if (rmid2
< rmid1
) rhi
= rhi
+ 1; /* carry? incr hi */
1071 r
->hi
= rhi
& LMASK
; /* mask result */
1076 void dp_neg (UDP
*r
)
1078 r
->lo
= NEG (r
->lo
);
1079 r
->hi
= (~r
->hi
+ (r
->lo
== 0)) & LMASK
;
1083 /* Support routines */
1085 void unpackf (uint32 hi
, UFP
*r
)
1087 r
->sign
= hi
& FPSIGN
; /* get sign */
1088 r
->exp
= FD_GETEXP (hi
); /* get exponent */
1089 if (r
->exp
== 0) { /* exp = 0? */
1090 if (r
->sign
) RSVD_OPND_FAULT
; /* if -, rsvd op */
1091 r
->frac
.hi
= r
->frac
.lo
= 0; /* else 0 */
1094 r
->frac
.hi
= WORDSWAP ((hi
& ~(FPSIGN
| FD_EXP
)) | FD_HB
);
1096 dp_lsh (&r
->frac
, FD_GUARD
);
1100 void unpackd (uint32 hi
, uint32 lo
, UFP
*r
)
1102 r
->sign
= hi
& FPSIGN
; /* get sign */
1103 r
->exp
= FD_GETEXP (hi
); /* get exponent */
1104 if (r
->exp
== 0) { /* exp = 0? */
1105 if (r
->sign
) RSVD_OPND_FAULT
; /* if -, rsvd op */
1106 r
->frac
.hi
= r
->frac
.lo
= 0; /* else 0 */
1109 r
->frac
.hi
= WORDSWAP ((hi
& ~(FPSIGN
| FD_EXP
)) | FD_HB
);
1110 r
->frac
.lo
= WORDSWAP (lo
);
1111 dp_lsh (&r
->frac
, FD_GUARD
);
1115 void unpackg (uint32 hi
, uint32 lo
, UFP
*r
)
1117 r
->sign
= hi
& FPSIGN
; /* get sign */
1118 r
->exp
= G_GETEXP (hi
); /* get exponent */
1119 if (r
->exp
== 0) { /* exp = 0? */
1120 if (r
->sign
) RSVD_OPND_FAULT
; /* if -, rsvd op */
1121 r
->frac
.hi
= r
->frac
.lo
= 0; /* else 0 */
1124 r
->frac
.hi
= WORDSWAP ((hi
& ~(FPSIGN
| G_EXP
)) | G_HB
);
1125 r
->frac
.lo
= WORDSWAP (lo
);
1126 dp_lsh (&r
->frac
, G_GUARD
);
1133 static uint32 normmask
[5] = {
1134 0xc0000000, 0xf0000000, 0xff000000, 0xffff0000, 0xffffffff
1136 static int32 normtab
[6] = { 1, 2, 4, 8, 16, 32};
1138 if ((r
->frac
.hi
== 0) && (r
->frac
.lo
== 0)) { /* if fraction = 0 */
1139 r
->sign
= r
->exp
= 0; /* result is 0 */
1142 while ((r
->frac
.hi
& UF_NM_H
) == 0) { /* normalized? */
1143 for (i
= 0; i
< 5; i
++) { /* find first 1 */
1144 if (r
->frac
.hi
& normmask
[i
]) break;
1146 dp_lsh (&r
->frac
, normtab
[i
]); /* shift frac */
1147 r
->exp
= r
->exp
- normtab
[i
]; /* decr exp */
1152 int32
rpackfd (UFP
*r
, int32
*rh
)
1154 static UDP f_round
= { UF_FRND_L
, UF_FRND_H
};
1155 static UDP d_round
= { UF_DRND_L
, UF_DRND_H
};
1157 if (rh
) *rh
= 0; /* assume 0 */
1158 if ((r
->frac
.hi
== 0) && (r
->frac
.lo
== 0)) return 0; /* result 0? */
1159 if (rh
) dp_add (&r
->frac
, &d_round
); /* round */
1160 else dp_add (&r
->frac
, &f_round
);
1161 if ((r
->frac
.hi
& UF_NM_H
) == 0) { /* carry out? */
1162 dp_rsh (&r
->frac
, 1); /* renormalize */
1163 r
->exp
= r
->exp
+ 1;
1165 if (r
->exp
> (int32
) FD_M_EXP
) FLT_OVFL_FAULT
; /* ovflo? fault */
1166 if (r
->exp
<= 0) { /* underflow? */
1167 if (PSL
& PSW_FU
) FLT_UNFL_FAULT
; /* fault if fu */
1168 return 0; /* else 0 */
1170 dp_rsh (&r
->frac
, FD_GUARD
); /* remove guard */
1171 if (rh
) *rh
= WORDSWAP (r
->frac
.lo
); /* get low */
1172 return r
->sign
| (r
->exp
<< FD_V_EXP
) |
1173 (WORDSWAP (r
->frac
.hi
) & ~(FD_HB
| FPSIGN
| FD_EXP
));
1176 int32
rpackg (UFP
*r
, int32
*rh
)
1178 static UDP g_round
= { UF_GRND_L
, UF_GRND_H
};
1180 *rh
= 0; /* assume 0 */
1181 if ((r
->frac
.hi
== 0) && (r
->frac
.lo
== 0)) return 0; /* result 0? */
1182 dp_add (&r
->frac
, &g_round
); /* round */
1183 if ((r
->frac
.hi
& UF_NM_H
) == 0) { /* carry out? */
1184 dp_rsh (&r
->frac
, 1); /* renormalize */
1185 r
->exp
= r
->exp
+ 1;
1187 if (r
->exp
> (int32
) G_M_EXP
) FLT_OVFL_FAULT
; /* ovflo? fault */
1188 if (r
->exp
<= 0) { /* underflow? */
1189 if (PSL
& PSW_FU
) FLT_UNFL_FAULT
; /* fault if fu */
1190 return 0; /* else 0 */
1192 dp_rsh (&r
->frac
, G_GUARD
); /* remove guard */
1193 *rh
= WORDSWAP (r
->frac
.lo
); /* get low */
1194 return r
->sign
| (r
->exp
<< G_V_EXP
) |
1195 (WORDSWAP (r
->frac
.hi
) & ~(G_HB
| FPSIGN
| G_EXP
));
1200 /* Floating point instructions */
1202 /* Move/test/move negated floating
1204 Note that only the high 32b is processed.
1205 If the high 32b is not zero, it is unchanged.
1208 int32
op_movfd (int32 val
)
1210 if (val
& FD_EXP
) return val
;
1211 if (val
& FPSIGN
) RSVD_OPND_FAULT
;
1215 int32
op_mnegfd (int32 val
)
1217 if (val
& FD_EXP
) return (val
^ FPSIGN
);
1218 if (val
& FPSIGN
) RSVD_OPND_FAULT
;
1222 int32
op_movg (int32 val
)
1224 if (val
& G_EXP
) return val
;
1225 if (val
& FPSIGN
) RSVD_OPND_FAULT
;
1229 int32
op_mnegg (int32 val
)
1231 if (val
& G_EXP
) return (val
^ FPSIGN
);
1232 if (val
& FPSIGN
) RSVD_OPND_FAULT
;
1236 /* Floating to floating convert - F to D is essentially done with MOVFD */
1238 int32
op_cvtdf (int32
*opnd
)
1242 unpackd (opnd
[0], opnd
[1], &a
);
1243 return rpackfd (&a
, NULL
);
1246 int32
op_cvtfg (int32
*opnd
, int32
*rh
)
1250 unpackf (opnd
[0], &a
);
1251 a
.exp
= a
.exp
- FD_BIAS
+ G_BIAS
;
1252 return rpackg (&a
, rh
);
1255 int32
op_cvtgf (int32
*opnd
)
1259 unpackg (opnd
[0], opnd
[1], &a
);
1260 a
.exp
= a
.exp
- G_BIAS
+ FD_BIAS
;
1261 return rpackfd (&a
, NULL
);
1264 /* Floating add and subtract */
1266 int32
op_addf (int32
*opnd
, t_bool sub
)
1270 unpackf (opnd
[0], &a
); /* F format */
1271 unpackf (opnd
[1], &b
);
1272 if (sub
) a
.sign
= a
.sign
^ FPSIGN
; /* sub? -s1 */
1273 vax_fadd (&a
, &b
); /* add fractions */
1274 return rpackfd (&a
, NULL
);
1277 int32
op_addd (int32
*opnd
, int32
*rh
, t_bool sub
)
1281 unpackd (opnd
[0], opnd
[1], &a
);
1282 unpackd (opnd
[2], opnd
[3], &b
);
1283 if (sub
) a
.sign
= a
.sign
^ FPSIGN
; /* sub? -s1 */
1284 vax_fadd (&a
, &b
); /* add fractions */
1285 return rpackfd (&a
, rh
);
1288 int32
op_addg (int32
*opnd
, int32
*rh
, t_bool sub
)
1292 unpackg (opnd
[0], opnd
[1], &a
);
1293 unpackg (opnd
[2], opnd
[3], &b
);
1294 if (sub
) a
.sign
= a
.sign
^ FPSIGN
; /* sub? -s1 */
1295 vax_fadd (&a
, &b
); /* add fractions */
1296 return rpackg (&a
, rh
); /* round and pack */
1299 /* Floating multiply */
1301 int32
op_mulf (int32
*opnd
)
1305 unpackf (opnd
[0], &a
); /* F format */
1306 unpackf (opnd
[1], &b
);
1307 vax_fmul (&a
, &b
, 0, FD_BIAS
, 0, 0); /* do multiply */
1308 return rpackfd (&a
, NULL
); /* round and pack */
1311 int32
op_muld (int32
*opnd
, int32
*rh
)
1315 unpackd (opnd
[0], opnd
[1], &a
); /* D format */
1316 unpackd (opnd
[2], opnd
[3], &b
);
1317 vax_fmul (&a
, &b
, 1, FD_BIAS
, 0, 0); /* do multiply */
1318 return rpackfd (&a
, rh
); /* round and pack */
1321 int32
op_mulg (int32
*opnd
, int32
*rh
)
1325 unpackg (opnd
[0], opnd
[1], &a
); /* G format */
1326 unpackg (opnd
[2], opnd
[3], &b
);
1327 vax_fmul (&a
, &b
, 1, G_BIAS
, 0, 0); /* do multiply */
1328 return rpackg (&a
, rh
); /* round and pack */
1331 /* Floating divide */
1333 int32
op_divf (int32
*opnd
)
1337 unpackf (opnd
[0], &a
); /* F format */
1338 unpackf (opnd
[1], &b
);
1339 vax_fdiv (&a
, &b
, 26, FD_BIAS
); /* do divide */
1340 return rpackfd (&b
, NULL
); /* round and pack */
1343 int32
op_divd (int32
*opnd
, int32
*rh
)
1347 unpackd (opnd
[0], opnd
[1], &a
); /* D format */
1348 unpackd (opnd
[2], opnd
[3], &b
);
1349 vax_fdiv (&a
, &b
, 58, FD_BIAS
); /* do divide */
1350 return rpackfd (&b
, rh
); /* round and pack */
1353 int32
op_divg (int32
*opnd
, int32
*rh
)
1357 unpackg (opnd
[0], opnd
[1], &a
); /* G format */
1358 unpackg (opnd
[2], opnd
[3], &b
);
1359 vax_fdiv (&a
, &b
, 55, G_BIAS
); /* do divide */
1360 return rpackg (&b
, rh
); /* round and pack */
1363 /* Polynomial evaluation
1364 The most mis-implemented instruction in the VAX (probably here too).
1365 POLY requires a precise combination of masking versus normalizing
1366 to achieve the desired answer. In particular, the multiply step
1367 is masked prior to normalization. In addition, negative small
1368 fractions must not be treated as 0 during denorm.
1371 void op_polyf (int32
*opnd
, int32 acc
)
1374 int32 deg
= opnd
[1];
1375 int32 ptr
= opnd
[2];
1378 if (deg
> 31) RSVD_OPND_FAULT
; /* degree > 31? fault */
1379 unpackf (opnd
[0], &a
); /* unpack arg */
1380 wd
= Read (ptr
, L_LONG
, RD
); /* get C0 */
1382 unpackf (wd
, &r
); /* unpack C0 */
1383 res
= rpackfd (&r
, NULL
); /* first result */
1384 for (i
= 0; i
< deg
; i
++) { /* loop */
1385 unpackf (res
, &r
); /* unpack result */
1386 vax_fmul (&r
, &a
, 0, FD_BIAS
, 1, LMASK
); /* r = r * arg, mask */
1387 wd
= Read (ptr
, L_LONG
, RD
); /* get Cnext */
1389 unpackf (wd
, &c
); /* unpack Cnext */
1390 vax_fadd (&r
, &c
); /* r = r + Cnext */
1391 res
= rpackfd (&r
, NULL
); /* round and pack */
1399 void op_polyd (int32
*opnd
, int32 acc
)
1402 int32 deg
= opnd
[2];
1403 int32 ptr
= opnd
[3];
1404 int32 i
, wd
, wd1
, res
, resh
;
1406 if (deg
> 31) RSVD_OPND_FAULT
; /* degree > 31? fault */
1407 unpackd (opnd
[0], opnd
[1], &a
); /* unpack arg */
1408 wd
= Read (ptr
, L_LONG
, RD
); /* get C0 */
1409 wd1
= Read (ptr
+ 4, L_LONG
, RD
);
1411 unpackd (wd
, wd1
, &r
); /* unpack C0 */
1412 res
= rpackfd (&r
, &resh
); /* first result */
1413 for (i
= 0; i
< deg
; i
++) { /* loop */
1414 unpackd (res
, resh
, &r
); /* unpack result */
1415 vax_fmul (&r
, &a
, 1, FD_BIAS
, 0, 1); /* r = r * arg, mask */
1416 wd
= Read (ptr
, L_LONG
, RD
); /* get Cnext */
1417 wd1
= Read (ptr
+ 4, L_LONG
, RD
);
1419 unpackd (wd
, wd1
, &c
); /* unpack Cnext */
1420 vax_fadd (&r
, &c
); /* r = r + Cnext */
1421 res
= rpackfd (&r
, &resh
); /* round and pack */
1432 void op_polyg (int32
*opnd
, int32 acc
)
1435 int32 deg
= opnd
[2];
1436 int32 ptr
= opnd
[3];
1437 int32 i
, wd
, wd1
, res
, resh
;
1439 if (deg
> 31) RSVD_OPND_FAULT
; /* degree > 31? fault */
1440 unpackg (opnd
[0], opnd
[1], &a
); /* unpack arg */
1441 wd
= Read (ptr
, L_LONG
, RD
); /* get C0 */
1442 wd1
= Read (ptr
+ 4, L_LONG
, RD
);
1444 unpackg (wd
, wd1
, &r
); /* unpack C0 */
1445 res
= rpackg (&r
, &resh
); /* first result */
1446 for (i
= 0; i
< deg
; i
++) { /* loop */
1447 unpackg (res
, resh
, &r
); /* unpack result */
1448 vax_fmul (&r
, &a
, 1, G_BIAS
, 0, 1); /* r = r * arg */
1449 wd
= Read (ptr
, L_LONG
, RD
); /* get Cnext */
1450 wd1
= Read (ptr
+ 4, L_LONG
, RD
);
1452 unpackg (wd
, wd1
, &c
); /* unpack Cnext */
1453 vax_fadd (&r
, &c
); /* r = r + Cnext */
1454 res
= rpackg (&r
, &resh
); /* round and pack */