|
|
1.1 ! root 1: # double sqrt(arg):revised July 18,1980 ! 2: # double arg ! 3: # if(arg<0.0) { _errno=EDOM; return(-sqrt(-arg)) } ! 4: # W. Kahan's magic sqrt ! 5: # coded by Heidi Stettner ! 6: ! 7: .set EDOM,98 ! 8: .text ! 9: .align 1 ! 10: .globl _sqrt ! 11: .globl dsqrt_r5 ! 12: .globl _errno ! 13: _sqrt: ! 14: .word 0x003c # save r5,r4,r3,r2 ! 15: bispsw $0xe0 ! 16: movd 4(ap),r0 ! 17: bsbb dsqrt_r5 ! 18: ret ! 19: dsqrt_r5: ! 20: movd r0,r4 ! 21: jleq nonpos #argument is not positive ! 22: movzwl r4,r2 ! 23: ashl $-1,r2,r0 ! 24: addw2 $0x203c,r0 #r0 has magic initial appx ! 25: ! 26: # Do two steps of Heron's rule ! 27: ! 28: divf3 r0,r4,r2 ! 29: addf2 r2,r0 ! 30: subw2 $0x80,r0 ! 31: ! 32: divf3 r0,r4,r2 ! 33: addf2 r2,r0 ! 34: subw2 $0x80,r0 ! 35: ! 36: ! 37: # Scale argument and approximation to prevent over/underflow ! 38: # NOTE: The following four steps would not be necessary if underflow ! 39: # were gentle. ! 40: ! 41: bicw3 $0xffff807f,r4,r1 ! 42: subw2 $0x4080,r1 # r1 contains scaling factor ! 43: subw2 r1,r4 ! 44: movl r0,r2 ! 45: subw2 r1,r2 ! 46: ! 47: # Cubic step ! 48: ! 49: clrl r1 ! 50: clrl r3 ! 51: muld2 r0,r2 ! 52: subd2 r2,r4 ! 53: addw2 $0x100,r2 ! 54: addd2 r4,r2 ! 55: muld2 r0,r4 ! 56: divd2 r2,r4 ! 57: addw2 $0x80,r4 ! 58: addd2 r4,r0 ! 59: rsb ! 60: nonpos: ! 61: jneq negarg ! 62: clrd r0 #argument is zero ! 63: rsb ! 64: negarg: ! 65: movl $EDOM,_errno ! 66: mnegd r4,-(sp) ! 67: calls $2,_sqrt ! 68: mnegd r0,r0 # returns -sqrt(-arg) ! 69: ret
This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.