Annotation of researchv10dc/libnm/sqrt.s, revision 1.1.1.1

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:        movd    4(ap),r0
                     16:        bsbb    dsqrt_r5
                     17:        ret
                     18: dsqrt_r5:
                     19:        movd    r0,r4
                     20:        jleq    nonpos          #argument is not positive
                     21:        movzwl  r4,r2
                     22:        ashl    $-1,r2,r0
                     23:        addw2   $0x203c,r0      #r0 has magic initial appx
                     24: 
                     25: # Do two steps of Heron's rule
                     26: 
                     27:        divf3   r0,r4,r2        
                     28:        addf2   r2,r0
                     29:        subw2   $0x80,r0
                     30: 
                     31:        divf3   r0,r4,r2
                     32:        addf2   r2,r0
                     33:        subw2   $0x80,r0
                     34: 
                     35: 
                     36: # Scale argument and approximation to prevent over/underflow
                     37: # NOTE: The following four steps would not be necessary if underflow
                     38: #       were gentle.
                     39: 
                     40:        bicw3   $0xffff807f,r4,r1
                     41:        subw2   $0x4080,r1              # r1 contains scaling factor
                     42:        subw2   r1,r4
                     43:        movl    r0,r2
                     44:        subw2   r1,r2
                     45: 
                     46: # Cubic step
                     47: 
                     48:        clrl    r1
                     49:        clrl    r3
                     50:        muld2   r0,r2
                     51:        subd2   r2,r4
                     52:        addw2   $0x100,r2
                     53:        addd2   r4,r2
                     54:        muld2   r0,r4
                     55:        divd2   r2,r4
                     56:        addw2   $0x80,r4
                     57:        addd2   r4,r0
                     58:        rsb
                     59: nonpos:
                     60:        jneq    negarg
                     61:        clrd    r0              #argument is zero
                     62:        rsb
                     63: negarg:
                     64:        movl    $EDOM,_errno
                     65:        mnegd   r4,-(sp)
                     66:        calls   $2,_sqrt
                     67:        mnegd   r0,r0           # returns -sqrt(-arg)
                     68:        ret

unix.superglobalmegacorp.com

This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.