1*34128Sbostic/* 231840Szliu * Copyright (c) 1987 Regents of the University of California. 3*34128Sbostic * All rights reserved. 4*34128Sbostic * 5*34128Sbostic * Redistribution and use in source and binary forms are permitted 6*34128Sbostic * provided that this notice is preserved and that due credit is given 7*34128Sbostic * to the University of California at Berkeley. The name of the University 8*34128Sbostic * may not be used to endorse or promote products derived from this 9*34128Sbostic * software without specific prior written permission. This software 10*34128Sbostic * is provided ``as is'' without express or implied warranty. 11*34128Sbostic * 12*34128Sbostic * All recipients should regard themselves as participants in an ongoing 13*34128Sbostic * research project and hence should feel obligated to report their 14*34128Sbostic * experiences (good or bad) with these elementary function codes, using 15*34128Sbostic * the sendbug(8) program, to the authors. 1631840Szliu */ 1731840Szliu .data 1831840Szliu .align 2 1931840Szliu_sccsid: 20*34128Sbostic .asciz "@(#)support.s 5.3 (ucb.elefunt) 04/29/88" 2131840Szliu/* 2231840Szliu * copysign(x,y), 2331840Szliu * logb(x), 2431840Szliu * scalb(x,N), 2531840Szliu * finite(x), 2631840Szliu * drem(x,y), 2731840Szliu * Coded in vax assembly language by K. C. Ng 4/9/85. 2831840Szliu * Re-coded in tahoe assembly language by Z. Alex Liu 7/13/87. 2931840Szliu */ 3031840Szliu/* 3131840Szliu * double copysign(x,y) 3231840Szliu * double x,y; 3331840Szliu */ 3431840Szliu .globl _copysign 3531840Szliu .text 3631842Szliu .align 2 3731840Szliu_copysign: 3831840Szliu .word 0x0004 # save r2 3931840Szliu movl 8(fp),r1 4031840Szliu movl 4(fp),r0 # r0:r1 = x 4131840Szliu andl3 $0x7f800000,r0,r2 # r2 = biased exponent of x 4231840Szliu beql 1f # if 0 or reserved op then return x 4331840Szliu andl3 $0x80000000,12(fp),r2 # r2 = sign bit of y at bit-31 4431840Szliu andl2 $0x7fffffff,r0 # replace x by |x| 4531840Szliu orl2 r2,r0 # copy the sign bit of y to x 4631840Szliu1: ret 4731840Szliu/* 4831840Szliu * double logb(x) 4931840Szliu * double x; 5031840Szliu */ 5131840Szliu .globl _logb 5231840Szliu .text 5331842Szliu .align 2 5431840Szliu_logb: 5531842Szliu .word 0x0000 # save nothing 5631840Szliu andl3 $0x7f800000,4(fp),r0 # r0[b23:b30] = biased exponent of x 5731840Szliu beql 1f 5831840Szliu shrl $23,r0,r0 # r0[b0:b7] = biased exponent of x 5931840Szliu subl2 $129,r0 # r0 = unbiased exponent of x 6031840Szliu cvld r0 # acc = unbiased exponent of x (double) 6131840Szliu std r0 # r0 = unbiased exponent of x (double) 6231840Szliu ret 6331840Szliu1: movl 8(fp),r1 # 8(fp) must be moved first 6431840Szliu movl 4(fp),r0 # r0:r1 = x (zero or reserved op) 6531840Szliu blss 2f # simply return if reserved op 6631840Szliu movl $0xfe000000,r1 6731840Szliu movl $0xcfffffff,r0 # -2147483647.0 6831840Szliu2: ret 6931840Szliu/* 7031840Szliu * long finite(x) 7131840Szliu * double x; 7231840Szliu */ 7331840Szliu .globl _finite 7431840Szliu .text 7531842Szliu .align 2 7631840Szliu_finite: 7731842Szliu .word 0x0000 # save nothing 7831840Szliu andl3 $0xff800000,4(fp),r0 # r0 = sign of x & its biased exponent 7931840Szliu cmpl r0,$0x80000000 # is x a reserved op? 8031840Szliu beql 1f # if so, return FALSE (0) 8131840Szliu movl $1,r0 # else return TRUE (1) 8231840Szliu ret 8331840Szliu1: clrl r0 8431840Szliu ret 8531840Szliu/* 8631840Szliu * double scalb(x,N) 8731840Szliu * double x; int N; 8831840Szliu */ 8931840Szliu .globl _scalb 9031840Szliu .set ERANGE,34 9131840Szliu .text 9231842Szliu .align 2 9331840Szliu_scalb: 9431842Szliu .word 0x000c # save r2-r3 9531840Szliu movl 8(fp),r1 9631840Szliu movl 4(fp),r0 # r0:r1 = x (-128 <= Ex <= 126) 9731840Szliu andl3 $0x7f800000,r0,r3 # r3[b23:b30] = biased exponent of x 9831840Szliu beql 1f # is x a 0 or a reserved operand? 9931840Szliu movl 12(fp),r2 # r2 = N 10031840Szliu cmpl r2,$0xff # if N >= 255 10131840Szliu bgeq 2f # then the result must overflow 10231840Szliu cmpl r2,$-0xff # if N <= -255 10331840Szliu bleq 3f # then the result must underflow 10431840Szliu shrl $23,r3,r3 # r3[b0:b7] = biased exponent of x 10531840Szliu addl2 r2,r3 # r3 = biased exponent of the result 10631840Szliu bleq 3f # if <= 0 then the result underflows 10731840Szliu cmpl r3,$0x100 # if >= 256 then the result overflows 10831840Szliu bgeq 2f 10931840Szliu shll $23,r3,r3 # r3[b23:b30] = biased exponent of res. 11031840Szliu andl2 $0x807fffff,r0 11131840Szliu orl2 r3,r0 # r0:r1 = x*2^N 11231840Szliu1: ret 11331840Szliu2: pushl $ERANGE # if the result would overflow 11431840Szliu callf $8,_infnan # and _infnan returns 11531840Szliu andl3 $0x80000000,4(fp),r2 # get the sign of input arg 11631840Szliu orl2 r2,r0 # re-attach the sign to r0:r1 11731840Szliu ret 11831840Szliu3: clrl r1 # if the result would underflow 11931840Szliu clrl r0 # then return 0 12031840Szliu ret 12131840Szliu/* 12231840Szliu * double drem(x,y) 12331840Szliu * double x,y; 12431840Szliu * Returns x-n*y where n=[x/y] rounded (to even in the half way case). 12531840Szliu */ 12631840Szliu .globl _drem 12731840Szliu .set EDOM,33 12831840Szliu .text 12931842Szliu .align 2 13031840Szliu_drem: 13131840Szliu .word 0x1ffc # save r2-r12 13231840Szliu movl 16(fp),r3 13331840Szliu movl 12(fp),r2 # r2:r3 = y 13431840Szliu movl 8(fp),r1 13531840Szliu movl 4(fp),r0 # r0:r1 = x 13631840Szliu andl3 $0xff800000,r0,r4 13731840Szliu cmpl r4,$0x80000000 # is x a reserved operand? 13831840Szliu beql 1f # if yes then propagate x and return 13931840Szliu andl3 $0xff800000,r2,r4 14031840Szliu cmpl r4,$0x80000000 # is y a reserved operand? 14131840Szliu bneq 2f 14231840Szliu movl r3,r1 14331840Szliu movl r2,r0 # if yes then propagate y and return 14431840Szliu1: ret 14531840Szliu 14631840Szliu2: tstl r4 # is y a 0? 14731840Szliu bneq 3f 14831840Szliu pushl $EDOM # if so then generate reserved op fault 14931840Szliu callf $8,_infnan 15031840Szliu ret 15131840Szliu 15231840Szliu3: andl2 $0x7fffffff,r2 # r2:r3 = y <- |y| 15331840Szliu clrl r12 # r12 = nx := 0 15431840Szliu cmpl r2,$0x1c800000 # Ey ? 57 15531840Szliu bgtr 4f # if Ey > 57 goto 4 15631840Szliu addl2 $0x1c800000,r2 # scale up y by 2**57 15731840Szliu movl $0x1c800000,r12 # r12[b23:b30] = nx = 57 15831840Szliu4: pushl r12 # pushed onto stack: nf := nx 15931840Szliu andl3 $0x80000000,r0,-(sp) # pushed onto stack: sign of x 16031840Szliu andl2 $0x7fffffff,r0 # r0:r1 = x <- |x| 16131840Szliu movl r3,r11 # r10:r11 = y1 = y w/ last 27 bits 0 16231840Szliu andl3 $0xf8000000,r10,r11 # clear last 27 bits of y1 16331840Szliu 16431840SzliuLoop: cmpd2 r0,r2 # x ? y 16531840Szliu bleq 6f # if x <= y goto 6 16631840Szliu /* # begin argument reduction */ 16731840Szliu movl r3,r5 16831840Szliu movl r2,r4 # r4:r5 = t = y 16931840Szliu movl r11,r7 17031840Szliu movl r10,r6 # r6:r7 = t1 = y1 17131840Szliu andl3 $0x7f800000,r0,r8 # r8[b23:b30] = Ex:biased exponent of x 17231840Szliu andl3 $0x7f800000,r2,r9 # r9[b23:b30] = Ey:biased exponent of y 17331840Szliu subl2 r9,r8 # r8[b23:b30] = Ex-Ey 17431840Szliu subl2 $0x0c800000,r8 # r8[b23:b30] = k = Ex-Ey-25 17531840Szliu blss 5f # if k < 0 goto 5 17631840Szliu addl2 r8,r4 # t += k 17731840Szliu addl2 r8,r6 # t1 += k, scale up t and t1 17831840Szliu5: ldd r0 # acc = x 17931840Szliu divd r4 # acc = x/t 18031840Szliu cvdl r8 # r8 = n = [x/t] truncated 18131840Szliu cvld r8 # acc = dble(n) 18231840Szliu std r8 # r8:r9 = dble(n) 18331840Szliu ldd r4 # acc = t 18431840Szliu subd r6 # acc = t-t1 18531840Szliu muld r8 # acc = n*(t-t1) 18631840Szliu std r4 # r4:r5 = n*(t-t1) 18731840Szliu ldd r6 # acc = t1 18831840Szliu muld r8 # acc = n*t1 18931840Szliu subd r0 # acc = n*t1-x 19031840Szliu negd # acc = x-n*t1 19131840Szliu subd r4 # acc = (x-n*t1)-n*(t-t1) 19231840Szliu std r0 # r0:r1 = (x-n*t1)-n*(t-t1) 19331840Szliu brb Loop 19431840Szliu 19531840Szliu6: movl r12,r6 # r6 = nx 19631840Szliu beql 7f # if nx == 0 goto 7 19731840Szliu addl2 r6,r0 # x <- x*2**57:scale x up by nx 19831840Szliu clrl r12 # clear nx 19931840Szliu brb Loop 20031840Szliu 20131840Szliu7: movl r3,r5 20231840Szliu movl r2,r4 # r4:r5 = y 20331840Szliu subl2 $0x800000,r4 # r4:r5 = y/2 20431840Szliu cmpd2 r0,r4 # x ? y/2 20531840Szliu blss 9f # if x < y/2 goto 9 20631840Szliu bgtr 8f # if x > y/2 goto 8 20731840Szliu ldd r8 # acc = dble(n) 20831840Szliu cvdl r8 # r8 = ifix(dble(n)) 20931840Szliu bbc $0,r8,9f # if the last bit is zero, goto 9 21031840Szliu8: ldd r0 # acc = x 21131840Szliu subd r2 # acc = x-y 21231840Szliu std r0 # r0:r1 = x-y 21331840Szliu9: xorl2 (sp)+,r0 # x^sign (exclusive or) 21431840Szliu movl (sp)+,r6 # r6 = nf 21531840Szliu andl3 $0x7f800000,r0,r8 # r8 = biased exponent of x 21631840Szliu andl2 $0x807fffff,r0 # r0 = x w/ exponent zapped 21731840Szliu subl2 r6,r8 # r8 = Ex-nf 21831840Szliu bgtr 0f # if Ex-nf > 0 goto 0 21931840Szliu clrl r8 # clear r8 22031840Szliu clrl r0 22131840Szliu clrl r1 # x underflows to zero 22231840Szliu0: orl2 r8,r0 # put r8 into x's exponent field 22331840Szliu ret 224